{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "plot2d (generic function with 2 methods)"
      ]
     },
     "execution_count": 1,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "using JuLIP, JuLIP.ASE, JuLIP.Potentials \n",
    "using PyPlot\n",
    "\n",
    "function plot2d(at, ttl = nothing)\n",
    "    x, y, _ = xyz(at)\n",
    "    plot(x, y, \"b.\", markersize=20)\n",
    "    if ttl != nothing\n",
    "        title(ttl)\n",
    "    end \n",
    "    PyPlot.draw()\n",
    "    PyPlot.pause(0.0001)    \n",
    "end \n",
    "\n",
    "# NOTES \n",
    "#\n",
    "# 1. For now, this notebook is just to implement a temporary\n",
    "#    example for MD, Newtonian dynamics. This functionality\n",
    "#    should be moved into the JuLIP library, using more\n",
    "#    sophisticated integrators. \n",
    "#\n",
    "# 2. PyPlot is horrendously slow, should maybe switch to GR\n",
    "#    or Plotly for generating the animation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAhoAAAGgCAYAAADsAM6oAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAAPYQAAD2EBqD+naQAAIABJREFUeJzt3X9sXfV9//HXuTZJzPW1jW/zg8xO41xZpQVasYWgkmqbRASqIja6rtUqChGVKg0lcUImFFiVrRKDNJ3axQSWFqQxR2thk5J0HRRbGaXpUEkTyGjL2oFHAk3KCHPie68dBxf5fr5/HDnfOCTxvefez/ncez7Ph3RlentOnufcn29fn3tvYIwxAgAAsCDlegMAAEByMWgAAABrGDQAAIA1DBoAAMAaBg0AAGANgwYAALCGQQMAAFjDoAEAAKxh0AAAANYwaAAAAGsYNAAAgDXNcQdLpZLefvttZTIZBUEQdx4AAERgjNHY2JgWL16sVKr81yliHzTefvttdXd3x50FAAA1cOzYMXV1dZW9fOyDRiaTkRRuaFtbW9x5AAAQQbFYVHd399nn8XLFPmhM/7mkra2NQQMAgAZT6WEPHAwKAACsYdAAAADWMGgAAABrGDQAAIA1DBoAAMCa2N91AgC+M0Y6eVIaH5daW6VsVuLzC5FUvKIBXIQx0siI9Oab4U9jkt920fXpcs7npf5+qbdXmj9f6ukJf/b2hufn83b7Pt2u6MZ7X7okE7NCoWAkmUKhEHcaKMvoqDHbtxuTyxkT3lXDUy4Xnj86mry2i65vl/PgoDHptDFBEJ7O7U6fl06Hy9WaT7cruva6UZ+/GTSAc7h8MnDVdtH17XIeHDSmqcmYVGpm7/xTKhUul4TLmm7yugwaQJVcPxm4aLvo+nY5j46GD/SzNc9tp9O1+S3Up9sVXfvdqM/fgTHGxPmnmmKxqPb2dhUKBT6CHHUjn5e6uqQzZ6RSafblUymppUU6flzq6GjMtouuj5dzf790zz3hQ365gkDavl3q64ve9el2RTeebtTnbw4GBSQNDEgTE+XdYaVwuYkJadeuxm276Pp2ORsj7dgRbd2HH65sODmfT7cruvF1o+AVDXjPmPCo/yNHKv+tc9kyaXg4+lsTXbVddH28nEdGwneVRDUyEr71tVI+3a7oxteN+vzNoAHvuXoycNl20fXxcn7zzfAtrFEdPSotXVr5ej7drujG1+VPJ0BE4+PVrT821nhtF10fL+fW1uq6mUy09Xy6XdGNrxsVgwa85+rJwGXbRdfHyzmblXK5yl+mDoJwvc7OaF2fbld04+tGxaAB77l6MnDZdtH18XIOAmn9+mjr9vVFPybFp9sV3fi6UTFowHuungxctl10fbycJWnNGunyy8O3F5YjlQqXv/PO6E2fbld04+tGxcGggPz8fAc+RyO+9tCQtHp1+A6BS7VTqfBJ4Ac/kG6+ubqmT7cruvF0ORgUqEJHh7R7d/ggP9tvntNPBnv2VP8E5LLtouvj5SxJt9wiPfNM+EAfBB/8jXL6vJaW2gwZkl+3K7rxdSOp7gNJK8dHkKOelfu9AUNDyWm76Pp4ORsTfqx4f/+Fv/yqv9+YfL72TZ9uV3TtdvmuE6BGXDwZuG676Pp4OU8rlYwZGTHm6NHwZ6lkt+fT7YquvS7fdQLUmDHSqVPhe84zmfBI7bgOonLVdtH18XJ2xafbFd3ad/lkUAAAYA0HgwIAgLrDoAEAAKxh0AAAANYwaAAAAGsYNAAAgDUMGgAAwBoGDQAAYA2DBgAAsIZBAwAAWMOgAQAArGHQAAAA1jBoAAAAayoaNKamprRlyxb19PSopaVFuVxODzzwgGL+XjYAANAgmitZeNu2bdq5c6cGBgZ09dVX66WXXtJdd92l9vZ29fX12dpGAADQoCoaNH7yk5/oj//4j7V69WpJ0tKlS/Xkk0/q4MGDF11ncnJSk5OTZ/93sViMuKkAAKDRVPSnkxtvvFHPPfecXn/9dUnSz372M73wwgv69Kc/fdF1tm7dqvb29rOn7u7u6rYYqIAx0siI9Oab4c+4/srnquuyTZcu3cbtWmUqMDU1ZTZv3myCIDDNzc0mCALz0EMPXXKd9957zxQKhbOnY8eOGUmmUChUkgYqMjpqzPbtxuRyxoR31fCUy4Xnj44mq+uyTZcu3cbtVqJQKER6/q5o0HjyySdNV1eXefLJJ83Pf/5zs2vXLtPZ2Wn+8R//0fqGAuUaHDQmnTYmCMLTuXfa6fPS6XC5JHRdtunSpdu43UrFMmh0dXWZRx55ZMZ5DzzwgPnIRz5S9r/BoAGbBgeNaWoyJpWaeWc9/5RKhcvV6o7rquuyTZcu3cbtRhHLoNHZ2Wn+/u//fsZ5Dz30kOnt7S3732DQgC2jo+HUP9sd9tw7bjpd/UuSrro+7jNdunTd/Rkl6vN3RQeD3nrrrXrwwQf1zDPP6M0339TevXv1zW9+U5/5zGdqfegIULGBAWliQiqVylu+VAqX37WrMbsu23Tp0m3cbuwqmUqKxaLZsGGDWbJkiZk3b55ZtmyZ+cpXvmImJyetT0TApZRK4UFT5/99c7ZTEITrlUqN1fVxn+nSpVubx46ooj5/B8YYE+dgUywW1d7erkKhoLa2tjjTSLCREWn+/OrWz2Ybp+uyTZcu3cbtViPq8zffdYJEGB+vbv2xscbqumzTpUu3cbsuMGggEVpbq1s/k2msrss2Xbp0G7frAoMGEiGblXI5KQgqWy8IwvU6Oxur67JNly7dxu26wKCBRAgCaf36aOv29VV+Z3fdddmmS5du43Zd4GBQJEY+L3V1SWfOlPd2sVRKammRjh+XOjoar+uyTZcu3cbtRsXBoPBeR4e0e3c46admuWWnUuFye/ZUf4d11XXZpkuXbuN2Y2fhrbaXxOdowLZyvzdgaCgZXZdtunTpNm63UrF8BHktMGggDqOjxvT3X/ibEPv7jcnnk9V12aZLl27jdivBB3YBF2CMdOpU+J7zTCY8UjuOg6hcdV226dKl27jdckR9/mbQAAAAs+JgUAAAUHcYNAAAgDUMGgAAwBoGDQAAYA2DBgAAsIZBAwAAWMOgAQAArGHQAAAA1jBoAAAAaxg0AACANQwaAADAGgYNAABgDYMGAACwhkEDAABYw6ABAACsYdAAAADWNLveAGA2xkgnT0rj41Jrq5TNSkGQ7DZdunTpJgWvaKBu5fNSf7/U2yvNny/19IQ/e3vD8/P55LXp0qVLN3FMzAqFgpFkCoVC3Gk0kMFBY9JpY4IgPIW/I4Sn6fPS6XC5pLTp0qVLt55Fff5m0EDdGRw0pqnJmFRq5p31/FMqFS5XyzuuqzZdunTp1ruoz9+BMcbE+QpKsVhUe3u7CoWC2tra4kyjAeTzUleXdOaMVCrNvnwqJbW0SMePSx0djdmmS5cu3UYQ9fmbYzRQVwYGpImJ8u6wUrjcxIS0a1fjtunSpUs3yXhFA3XDmPCgqSNHwv8uVxBIy5ZJw8PRj+x21aZLly7dRhH1+ZtBA3VjZCQ8Qrua9bPZxmrTpUuXbqPgTydoeOPj1a0/NtZ4bbp06dJNOgYN1I3W1urWz2Qar02XLl26SceggbqRzUq5XOV/swyCcL3OzsZr06VLl27SMWigbgSBtH59tHX7+qo7qMpVmy5dunSTjoNBUVf4HA26dOnSrU8cDIpE6OiQdu8OJ/3ULLfOVCpcbs+e2txhXbXp0qVLN9Fq/Amls+IjyFGOcr83YGgoOW26dOnSrWd81wkSZ3TUmP5+Y3K5mXfaXC48P59PXpsuXbp06xXfdYLEMkY6dSp8z3kmEx6pHddBVK7adOnSpVtv+GRQAABgDQeDAgCAusOgAQAArGHQAAAA1jBoAAAAaxg0AACANQwaAADAGgYNAABgDYMGAACwhkEDAABYw6ABAACsYdAAAADWMGgAAABrGDQAAIA1DBoAAMAaBg0AAGBNs+sNgB+MkU6elMbHpdZWKZuVgoBuUrou23TpJqmbRLyiAavyeam/X+rtlebPl3p6wp+9veH5+TzdRu66bNOlm6RuopmYFQoFI8kUCoW404jZ4KAx6bQxQRCewt8RwtP0eel0uBzdxuu6bNOlm6Ruo4j6/M2gASsGB41pajImlZp5Zz3/lEqFy9Xqjks3nq7LNl26Seo2kqjP34ExxsT5CkqxWFR7e7sKhYLa2triTCMm+bzU1SWdOSOVSrMvn0pJLS3S8eNSRwfdeu+6bNOlm6Ruo4n6/M0xGqi5gQFpYqK8O6wULjcxIe3aRbcRui7bdOkmqesLXtFATRkTHjR15Ej43+UKAmnZMml4ONqR3XTj6bps06WbpG4jiu0Vjd/85jf64he/qGw2q5aWFl177bV66aWXKv1nkFAnT0pvvFHZHVYKl3/jDenUKbr13HXZpks3SV2fVDRojI6OauXKlbrsssv07LPP6pe//KW+8Y1v6IorrrC1fWgw4+PVrT82Rreeuy7bdOkmqeuTij6wa9u2beru7tYTTzxx9ryenp6abxQaV2trdetnMnTrueuyTZdukro+qegVje9///tavny5Pve5z2nBggW67rrr9Pjjj19yncnJSRWLxRknJFc2K+Vylf/NMgjC9To76dZz12WbLt0kdX1S0aBx5MgR7dy5U729vRoaGtLdd9+tvr4+DQwMXHSdrVu3qr29/eypu7u76o1G/QoCaf36aOv29UU/qIpuPF2Xbbp0k9T1SUXvOpkzZ46WL1+un/zkJ2fP6+vr06FDh/Tiiy9ecJ3JyUlNTk6e/d/FYlHd3d286yTBfHsvvG9dl226dJPUbTSxvOvkyiuv1Mc+9rEZ5330ox/Vr3/964uuM3fuXLW1tc04Idk6OqTdu8NJPzXLLSyVCpfbs6f6OyzdeLou23TpJqnrjUo+RvQLX/iC+dSnPjXjvI0bN5pPfvKTZf8bfAS5P8r93oChIbqN2HXZpks3Sd1GEct3nRw8eNA0NzebBx980AwPD5vvfOc75vLLLzf/9E//ZH1D0ZhGR43p7zcml5t5p83lwvPzebqN3HXZpks3Sd1GENt3nTz99NO6//77NTw8rJ6eHm3atElf/vKXy16fTwb1kzHhB9uMjYVvB+vsjOcgKrrxHazm2z7TpeubqM/ffAQ5AACYFV+qBgAA6g6DBgAAsIZBAwAAWMOgAQAArGHQAAAA1jBoAAAAaxg0AACANQwaAADAGgYNAABgDYMGAACwhkEDAABYw6ABAACsYdAAAADWMGgAAABrGDQAAIA1DBoAAMCaZtcbgMZhjHTypDQ+LrW2StmsFAR06TZumy5d2McrGphVPi/190u9vdL8+VJPT/iztzc8P5+nS7ex2nTpIkYmZoVCwUgyhUIh7jQiGBw0Jp02JgjCU/g7QniaPi+dDpejS7cR2nTp2rpNJ13U528GDVzU4KAxTU3GpFIz76znn1KpcLla3XHpJrvrsk2Xrq3btA+iPn8HxhgT5ysoxWJR7e3tKhQKamtrizONCuTzUleXdOaMVCrNvnwqJbW0SMePSx0ddOnWX5suXRtdn0R9/uYYDVzQwIA0MVHeHVYKl5uYkHbtoku3Ptt06droYna8ooEPMCY8aOrIkfC/yxUE0rJl0vBwtCO76Sa767JNl66Nrm+iPn8zaOADRkbCI7SrWT+bpUu3ftp06dro+oY/naBmxserW39sjC7d+mrTpWuji/IwaOADWlurWz+ToUu3vtp06droojwMGviAbFbK5Sr/m2UQhOt1dtKlW19tunRtdFEeBg18QBBI69dHW7evL/pBVXST3XXZpkvXRhfl4WBQXJBv74WnG0/XZZsuXRtdn3AwKGqqo0PavTuc9FOz3EpSqXC5PXuqv8PSTXbXZZsuXRtdlKHGn1A6Kz6CvLGU+70BQ0N06TZGmy5dW7fppOO7TmDN6Kgx/f3G5HIz77S5XHh+Pk+XbmO16dJF5fiuE1hnjHTqVPie80wmPFI7joOo6Ca767JNly7KxyeDAgAAazgYFAAA1B0GDQAAYA2DBgAAsIZBAwAAWMOgAQAArGHQAAAA1jBoAAAAaxg0AACANQwaAADAGgYNAABgDYMGAACwhkEDAABYw6ABAACsYdAAAADWMGgAAABrml1vAGCTMdLJk9L4uNTaKmWzUhC43iq7XO2zj5e1C1y/aDS8ouEZY6SREenNN8OfxiSzm89L/f1Sb680f77U0xP+7O0Nz8/n7fZdXM6u9tnHy9pFl+s33usXNWRiVigUjCRTKBTiTnttdNSY7duNyeWMCe+q4SmXC88fHU1Od3DQmHTamCAIT+d2p89Lp8Plas3V5exqn328rH26Tft4/eLioj5/M2h4wKcHqcFBY5qajEmlZvbOP6VS4XK1bru6nF3ss6+XtS+3aR+vX1wagwYuyKcHqdHR8MFntua57XS6Nr8ZubqcXe2zj5e1T7dpH69fzC7q8zfHaCRYPi999rPh3bJUuvSypVK43Gc/W/3fXF11BwakiYnZm+e2JyakXbuq67raX8ndPvt2Wft2m/bt+oVdDBoJ5tODlDHSjh3R1n344XD9qFxdzq722cfL2qfbtI/XL+wKjKnmZlG5YrGo9vZ2FQoFtbW1xZn2ijHhUeFHjlR2xw8CadkyaXg42lvXXHVHRsIj4aMaGQnfrlcpV/srudtn3y5r327Tvl2/KF/U529e0UiokyelN96o/LcLY8L1Tp1qrO74eLT1po2NRVvP1f5K7vbZt8vat9u0b9cv7GPQSCjfHqRaW6vrZjLR1nO1v5K7ffbtsvbtNu3b9Qv7GDQSyrcHqWxWyuUqf+k0CML1OjujdV3tr+Run327rH27Tft2/cI+Bo2E8u1BKgik9eujrdvXF/1vu672d/rfcLHPvl3Wvt2mfbt+YR+DRkL59iAlSWvWSJdfLqXKvFWnUuHyd94ZvelyfyU3++yqy216dly/qEtWPtXjEvjArvjwYT+zf9jP0FD1TZf7a4ybfXbV5TbN9WvzvoRL4wO78AEdHdLu3eGkP9tvJqlUuNyePeF6jdiVpFtukZ55RmppCf/d83/LmT6vpUX6wQ+km2+uvulyfyU3++yqy22a63eajfsSLLE0+FwUr2jEr9zvDajVb0Kuu8aEv+H091/4C5n6+43J52vfdLm/xrjZZ1ddbtNcvzbvS7gwvusEl+TTg9S5SiVjRkaMOXo0/Fkq2e253l9j4t9nV11u01y/iFfU528+GdQzxoQfbDM2Fr4drLMznoOoXHVd8W1/XeI2nWxczvXDySeDfu1rX1MQBNq4cWM1/wxiFATh28iWLg1/xnWHddV1xbf9dYnbdLJxOTe+yIPGoUOH9O1vf1sf//jHa7k9AAAgQSINGuPj47r99tv1+OOP64orrrjkspOTkyoWizNOAADAD5EGjbVr12r16tVatWrVrMtu3bpV7e3tZ0/d3d1RkgAAoAFVPGg89dRTOnz4sLZu3VrW8vfff78KhcLZ07FjxyreSAAA0JiaK1n42LFj2rBhg/bt26d58+aVtc7cuXM1d+7cSBsHAAAaW0Vvb/3e976nz3zmM2pqajp73tTUlIIgUCqV0uTk5Iz/70J4eysAAI0n6vN3Ra9o3HTTTfrFL34x47y77rpLV111lTZv3jzrkAEAAPxS0aCRyWR0zTXXzDgvnU4rm81+4HwAAAC+VA0AAFhT0SsaF/KjH/2oBpsBAACSiFc0AACANQwaAADAGgYNAABgDYMGAACwhkEDAABYU/W7ThA/Y6STJ6Xxcam1VcpmpSCgm6S2iy6XM92kdVEfeEWjgeTzUn+/1NsrzZ8v9fSEP3t7w/PzebqN3nbR5XKmm7Qu6oyJWaFQMJJMoVCIO93QBgeNSaeNCYLwFP6OEJ6mz0unw+XoNmbbRZfLmW7SurAn6vM3g0YDGBw0pqnJmFRq5p31/FMqFS5Xqzuub12XbRddLme6SevCrqjP3xV9e2st8O2tlcnnpa4u6cwZqVSafflUSmppkY4flzo66DZC20WXy5lu0rqwL+rzN8do1LmBAWliorw7rBQuNzEh7dpFt1HaLrpcznST1kX94hWNOmZMeNDUkSPhf5crCKRly6Th4WhHdvvWddl20eVyppu0LuIR9fmbQaOOjYyER2hXs342S7ee2y66XM50k9ZFPPjTSQKNj1e3/tgY3Xpvu+hyOdNNWhf1jUGjjrW2Vrd+JkO33tsuulzOdJPWRX1j0Khj2ayUy1X+N8sgCNfr7KRb720XXS5nuknror4xaNSxIJDWr4+2bl9f9IOqfOu6bLvocjnTTVoX9Y2DQeucb++F5/Md4ulyOdNNWhf2cTBoQnV0SLt3h5N+apZrK5UKl9uzp/o7rG9dl20XXS5nuknroo7V+BNKZ8VHkEdT7vcGDA3RbdS2iy6XM92kdWEP33XigdFRY/r7jcnlZt5pc7nw/HyebqO3XXS5nOkmrQs7+K4TjxgjnToVvuc8kwmP1I7jICrfui7bLrpcznST1kVt8cmgAADAGg4GBQAAdYdBAwAAWMOgAQAArGHQAAAA1jBoAAAAaxg0AACANQwaAADAGgYNAABgDYMGAACwhkEDAABYw6ABAACsYdAAAADWMGgAAABrGDQAAIA1DBoAAMCaZtcbANhkjHTypDQ+LrW2StmsFASutwq1xHUM1Dde0XDEGGlkRHrzzfCnMXRrKZ+X+vul3l5p/nyppyf82dsbnp/P2+27upxdtV00fb2OfesiAUzMCoWCkWQKhULc6bowOmrM9u3G5HLGhHfV8JTLheePjtKt1uCgMem0MUEQns7tTp+XTofL1Zqry9lV29X++ngd+9ZF/Yn6/M2gESNXD44+dQcHjWlqMiaVmtk7/5RKhcvVuu3qyc/VZe3qduXbdexbF/WJQaPOuXpw9Kk7Oho+6M3WPLedTtfmNzLXT35xt13tr4/XsW9d1C8GjTrm6sHRt+727R/8rWu2UxAY09/fmPvrqu1yf327jn3ror5Fff7mYNAYDAxIExNSqVTe8qVSuPyuXXTLZYy0Y0e0dR9+OFw/KleXs6u2q/318Tr2rYtkCoyp5u5XuWKxqPb2dhUKBbW1tcWZdsKY8Cj4I0cqe6ALAmnZMml4ONpb9XzrjoyE7ziIamQkfFtkpVztr6u2y/317Tr2rYv6F/X5m1c0LDt5Unrjjcp/mzImXO/UKbrlGB+Ptt60sbFo67naX1dtl/vr23XsWxfJxaBhmasHR9+6ra3VdTOZaOu52l9XbZf769t17FsXycWgYZmrB0ffutmslMtV/pJtEITrdXZG67raX1dtl/vr23XsWxfJxaBhmasHR9+6QSCtXx9t3b6+6H9TdrW/rtou99e369i3LpKLQcMyVw+OvnUlac0a6fLLpVSZt+pUKlz+zjujN13ur4u2y/2V/LqOfesiuXjXSQzyeamrSzpzpry3i6VSUkuLdPy41NFBtxJDQ9Lq1eGBaZdqp1LhA+IPfiDdfHN1TZf766Ltcn8lv65j37qob7zrpI51dEi7d4cPerP9Jjb94LhnT/V3WN+6knTLLdIzz4QPekHwwd+ups9raanNE5Dkdn9dtF3ur+TXdexbFwll4cPDLsnHTwadVu73BgwN0a3W6Gj4aZAX+iKo/n5j8vnaN13ur4u2y/01xq/r2Lcu6hMfQd4gXDw4+tidVioZMzJizNGj4c9SyW7P5f66aLu+fo3x5zr2rYv6E/X5m2M0HDEm/GCbsbHw7WCdnfEcROVb1xWX++ui7dv1K/l3X/LxOsZMUZ+/GTQAAMCsOBgUAADUHQYNAABgDYMGAACwhkEDAABYw6ABAACsYdAAAADWMGgAAABrGDQAAIA1DBoAAMCaigaNrVu36vrrr1cmk9GCBQt022236bXXXrO1bQAAoMFVNGjs379fa9eu1YEDB7Rv3z69//77uvnmm3X69Glb2wcAABpYVd918n//939asGCB9u/fr9///d8vax2+6wQAgMYT9fm7uZpooVCQJHV2dl50mcnJSU1OTp7938VisZokAABoIJEPBi2VStq4caNWrlypa6655qLLbd26Ve3t7WdP3d3dUZMAAKDBRP7Tyd13361nn31WL7zwgrq6ui663IVe0eju7uZPJwAANJBY/3Sybt06Pf300/rxj398ySFDkubOnau5c+dGydQ9Y6STJ6Xxcam1VcpmpSCgS7cxu67avu2vj134raI/nRhjtG7dOu3du1c//OEP1dPTY2u76lo+L/X3S7290vz5Uk9P+LO3Nzw/n6dLt3G6rtq+7a+PXUCSZCpw9913m/b2dvOjH/3I/O///u/Z08TERNn/RqFQMJJMoVCoJF03BgeNSaeNCYLwFP6OEJ6mz0unw+Xo0q33rqu2b/vrYxfJE/X5u6JBQ9IFT0888YT1Da0Hg4PGNDUZk0rNvLOef0qlwuVqdcelS9dG11Xbt/31sYtkivr8XdXnaETRqJ+jkc9LXV3SmTNSqTT78qmU1NIiHT8udXTQpVtfXVdt3/bXxy6SK+rzN991UqaBAWliorw7rBQuNzEh7dpFl279dV21fdtfH7vA+XhFowzGhAdNHTkS/ne5gkBatkwaHo52ZDdduja6rtq+7a+PXSRb1OdvBo0yjIyER2hXs342S5dufXRdtX3bXx+7SDb+dGLR+Hh164+N0aVbP11Xbd/218cucCEMGmVoba1u/UyGLt366bpq+7a/PnaBC2HQKEM2K+Vylf/NMgjC9S7xnXN06cbeddX2bX997AIXwqBRhiCQ1q+Ptm5fX/SDqujStdF11fZtf33sAhfCwaBl8u298HST3XXV9m1/fewiuTgY1LKODmn37nDST81yqaVS4XJ79lR/h6VL10bXVdu3/fWxC3xAjT+hdFaN/BHkxpT/vQFDQ3Tp1n/XVdu3/fWxi+SJ5btOaqHRBw1jjBkdNaa/35hcbuadNpcLz8/n6dJtnK6rtm/762MXycJ3nThgjHTqVPie80wmPFI7joOo6NJNUtu3/fWxi2Tgk0EBAIA1HAwKAADqDoMGAACwhkEDAABYw6ABAACsYdAAAADWMGgAAABrGDQAAIA1DBoAAMAaBg0AAGANgwYAALCGQQMegC2JAAAO9UlEQVQAAFjDoAEAAKxh0AAAANYwaAAAAGsYNAAAgDUMGgAAwJpm1xvgmjHSyZPS+LjU2ipls1IQ0KXbmF2XbbrJ7gJRefuKRj4v9fdLvb3S/PlST0/4s7c3PD+fp0u3cbou23ST3QWqZmJWKBSMJFMoFOJOnzU4aEw6bUwQhKfwd4TwNH1eOh0uR5duvXddtukmuwucK+rzt3eDxuCgMU1NxqRSM++s559SqXC5Wt1x6dK10XXZppvsLnC+qM/fgTHGxPkKSrFYVHt7uwqFgtra2uJMK5+XurqkM2ekUmn25VMpqaVFOn5c6uigS7e+ui7bdJPdBS4k6vO3V8doDAxIExPl3WGlcLmJCWnXLrp066/rsk032V2glrx5RcOY8KCpI0fC/y5XEEjLlknDw9GO7KZL10bXZZtusrvAxUR9/vZm0BgZCY/Qrmb9bJYu3froumzTTXYXuBj+dDKL8fHq1h8bo0u3frou23ST3QVqzZtBo7W1uvUzGbp066frsk032V2g1rwZNLJZKZer/G+WQRCu19lJl279dF226Sa7C9SaN4NGEEjr10dbt68v+kFVdOna6Lps0012F6g1bw4Glfx7LzzdZHddtukmuwtcCAeDlqGjQ9q9O5z0U7PseSoVLrdnT/V3WLp0bXRdtukmuwvUVI0/oXRWrj+C3JjyvzdgaIgu3frvumzTTXYXOBffdVKh0VFj+vuNyeVm3mlzufD8fJ4u3cbpumzTTXYXmMZ3nURkjHTqVPie80wmPFI7joOo6NJNWptusrsAnwwKAACs4WBQAABQdxg0AACANQwaAADAGgYNAABgDYMGAACwhkEDAABYw6ABAACsYdAAAADWMGgAAABrGDQAAIA1DBoAAMAaBg0AAGANgwYAALCGQQMAAFjDoAEAAKxpdr0BtWCMdPKkND4utbZK2awUBHTp0m2Ersu2b13AhYZ+RSOfl/r7pd5eaf58qacn/NnbG56fz9OlS7deuy7bvnUBp0zMCoWCkWQKhUJV/87goDHptDFBEJ7C3xHC0/R56XS4XC3RpUu3sdu+dYFaifr83ZCDxuCgMU1NxqRSM++s559SqXC5Wt1x6dKl29ht37pALUV9/g6MMSbOV1CKxaLa29tVKBTU1tZW8fr5vNTVJZ05I5VKsy+fSkktLdLx41JHR4QNpkuXbs26Ltu+dYFai/r8HekYjUcffVRLly7VvHnzdMMNN+jgwYNR/plIBgakiYny7rBSuNzEhLRrF126dF13XbZ96wJ1o9KXTp566ikzZ84c8w//8A/mv/7rv8yXv/xl09HRYU6cOGH1pRdjjCmVjMnlPvj3zdlOQRCuVypVnKRLl26Nuj7us8vLGqi12P50csMNN+j666/XI488IkkqlUrq7u7W+vXrdd999826fjV/OhkZCY/QjmpkJHwbWZT16NKlW13XZdu3LmBDLH86+e1vf6uXX35Zq1at+v//QCqlVatW6cUXX7zgOpOTkyoWizNOUY2PR15VkjQ2RpcuXVddl23fukA9qWjQGBkZ0dTUlBYuXDjj/IULF+qdd9654Dpbt25Ve3v72VN3d3fkjW1tjbyqJCmToUuXrquuy7ZvXaCeWP/Arvvvv1+FQuHs6dixY5H/rWxWyuUq/wS9IAjX6+ykS5euq67Ltm9doJ5UNGh86EMfUlNTk06cODHj/BMnTmjRokUXXGfu3Llqa2ubcYoqCKT166Ot29cX/SN+6dKlW33XZdu3LlBPIh0MumLFCu3YsUNSeDDokiVLtG7dOusHg0r+vReeLt0kdV22fesCtRbb52hs2rRJjz/+uAYGBvSrX/1Kd999t06fPq277rqr0n8qko4OaffucNJPzbL1qVS43J491d9h6dKlW5snPt/22eVlDdSFKO+l3bFjh1myZImZM2eOWbFihTlw4EDZ68b9XSdDQ1Vl6NKla6Hrsu1bF6gVr77rZNroqDH9/eEH25x7p83lwvPz+RpsMF26dK3xbZ9dXtZAtbz5rpMLMUY6dSp8z3kmEx6pHcdBVHTp0m3stm9doBpRn78TMWgAAAC7Yv1SNQAAgHIwaAAAAGsYNAAAgDUMGgAAwBoGDQAAYA2DBgAAsKY57uD0u2mLxWLcaQAAENH083aln4oR+6AxNjYmSeru7o47DQAAqjQ2Nqb29vayl4/9A7tKpZLefvttZTIZBTX8KLxisaju7m4dO3bMiw8CY3+Tjf1NPt/2mf1tfMYYjY2NafHixUrN9g2B54j9FY1UKqWuri5r/35bW1tirtRysL/Jxv4mn2/7zP42tkpeyZjGwaAAAMAaBg0AAGBN01e/+tWvut6IWmlqatIf/uEfqrk59r8IOcH+Jhv7m3y+7TP766fYDwYFAAD+4E8nAADAGgYNAABgDYMGAACwhkEDAABYw6ABAACsScSg8eijj2rp0qWaN2+ebrjhBh08eND1JlmzdetWXX/99cpkMlqwYIFuu+02vfbaa643KzZf+9rXFASBNm7c6HpTrPnNb36jL37xi8pms2ppadG1116rl156yfVmWTE1NaUtW7aop6dHLS0tyuVyeuCBByr+0qZ69eMf/1i33nqrFi9erCAI9L3vfW/G/2+M0V/91V/pyiuvVEtLi1atWqXh4WFHW1u9S+3v+++/r82bN+vaa69VOp3W4sWLdeedd+rtt992uMXVm+06Ptef//mfKwgCbd++PcYtdK/hB41//ud/1qZNm/TXf/3XOnz4sD7xiU/olltu0bvvvut606zYv3+/1q5dqwMHDmjfvn16//33dfPNN+v06dOuN826Q4cO6dvf/rY+/vGPu94Ua0ZHR7Vy5UpddtllevbZZ/XLX/5S3/jGN3TFFVe43jQrtm3bpp07d+qRRx7Rr371K23btk1f//rXtWPHDtebVhOnT5/WJz7xCT366KMX/P+//vWv6+GHH9a3vvUt/fSnP1U6ndYtt9yi9957L+YtrY1L7e/ExIQOHz6sLVu26PDhw9qzZ49ee+01/dEf/ZGDLa2d2a7jaXv37tWBAwe0ePHimLasjpgGt2LFCrN27dqz/3tqasosXrzYbN261eFWxefdd981ksz+/ftdb4pVY2Njpre31+zbt8/8wR/8gdmwYYPrTbJi8+bN5lOf+pTrzYjN6tWrzZe+9KUZ5/3Jn/yJuf322x1tkT2SzN69e8/+71KpZBYtWmT+9m//9ux5+XzezJ071zz55JMuNrGmzt/fCzl48KCRZN56662Ytsqui+3z8ePHze/8zu+YV1991Xz4wx82f/d3f+dg69xp6Fc0fvvb3+rll1/WqlWrzp6XSqW0atUqvfjiiw63LD6FQkGS1NnZ6XhL7Fq7dq1Wr14947pOou9///tavny5Pve5z2nBggW67rrr9Pjjj7veLGtuvPFGPffcc3r99dclST/72c/0wgsv6NOf/rTjLbPv6NGjeuedd2bcptvb23XDDTd49fgVBIE6Ojpcb4o1pVJJd9xxh+69915dffXVrjfHiYb+XNSRkRFNTU1p4cKFM85fuHCh/vu//9vRVsWnVCpp48aNWrlypa655hrXm2PNU089pcOHD+vQoUOuN8W6I0eOaOfOndq0aZP+8i//UocOHVJfX5/mzJmjNWvWuN68mrvvvvtULBZ11VVXqampSVNTU3rwwQd1++23u94069555x1JuuDj1/T/l2TvvfeeNm/erC984QuJ+nbT823btk3Nzc3q6+tzvSnONPSg4bu1a9fq1Vdf1QsvvOB6U6w5duyYNmzYoH379mnevHmuN8e6Uqmk5cuX66GHHpIkXXfddXr11Vf1rW99K5GDxr/8y7/oO9/5jr773e/q6quv1iuvvKKNGzdq8eLFidxfhN5//319/vOflzFGO3fudL051rz88svq7+/X4cOHFQSB681xpqH/dPKhD31ITU1NOnHixIzzT5w4oUWLFjnaqnisW7dOTz/9tJ5//nl1dXW53hxrXn75Zb377rv63d/9XTU3N6u5uVn79+/Xww8/rObmZk1NTbnexJq68sor9bGPfWzGeR/96Ef161//2tEW2XXvvffqvvvu05/92Z/p2muv1R133KF77rlHW7dudb1p1k0/Rvn2+DU9ZLz11lvat29fol/N+I//+A+9++67WrJkydnHr7feekt/8Rd/oaVLl7revNg09KAxZ84c/d7v/Z6ee+65s+eVSiU999xz+uQnP+lwy+wxxmjdunXau3evfvjDH6qnp8f1Jll100036Re/+IVeeeWVs6fly5fr9ttv1yuvvKKmpibXm1hTK1eu/MDblV9//XV9+MMfdrRFdk1MTCiVmvkw1NTUpFKp5GiL4tPT06NFixbNePwqFov66U9/mtjHr+khY3h4WP/+7/+ubDbrepOsuuOOO/Tzn/98xuPX4sWLde+992poaMj15sWm4f90smnTJq1Zs0bLly/XihUrtH37dp0+fVp33XWX602zYu3atfrud7+rf/3Xf1Umkzn7t9z29na1tLQ43rray2QyHzj+JJ1OK5vNJvK4lHvuuUc33nijHnroIX3+85/XwYMH9dhjj+mxxx5zvWlW3HrrrXrwwQe1ZMkSXX311frP//xPffOb39SXvvQl15tWE+Pj4/qf//mfs//76NGjeuWVV9TZ2aklS5Zo48aN+pu/+Rv19vaqp6dHW7Zs0eLFi3Xbbbc53OroLrW/V155pf70T/9Uhw8f1tNPP62pqamzj1+dnZ2aM2eOq82uymzX8fnD1GWXXaZFixbpIx/5SNyb6o7rt73Uwo4dO8ySJUvMnDlzzIoVK8yBAwdcb5I1ki54euKJJ1xvWmyS/PZWY4z5t3/7N3PNNdeYuXPnmquuuso89thjrjfJmmKxaDZs2GCWLFli5s2bZ5YtW2a+8pWvmMnJSdebVhPPP//8Be+va9asMcaEb3HdsmWLWbhwoZk7d6656aabzGuvveZ2o6twqf09evToRR+/nn/+edebHtls1/H5fHx7a2BMQj6CDwAA1J2GPkYDAADUNwYNAABgDYMGAACwhkEDAABYw6ABAACsYdAAAADWMGgAAABrGDQAAIA1DBoAAMAaBg0AAGANgwYAALDm/wF1oB8P/7ojngAAAABJRU5ErkJggg==",
      "text/plain": [
       "PyPlot.Figure(PyObject <matplotlib.figure.Figure object at 0x1351517d0>)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# problem parameters \n",
    "r0 = 1.0  # atomic spacing \n",
    "N = 11    # number of atoms in each direction \n",
    "\n",
    "# generate positions  (for real crystals this is automated)\n",
    "t = linspace(0.0, (N-1)*r0, N)\n",
    "o = ones(N)\n",
    "x, y = t * o', o * t'\n",
    "X = [ [1.0 0.5; 0.0 √3/2] * [x[:]'; y[:]']; zeros(N^2)' ]\n",
    "X = (X + 0.01*rand(size(X))) |> vecs\n",
    "# generate the atoms object  (H is just a place-holder)\n",
    "at = Atoms(\"H$(length(X))\", X) \n",
    "set_cell!(at, diagm([2.0*N, 2.0*N, 1.0]))\n",
    "set_pbc!(at, (true, true, true))\n",
    "set_constraint!(at, FixedCell(at))\n",
    "plot2d(at)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "# # generate a pair potential \n",
    "# lj = lennardjones() \n",
    "# set_calculator!(at, lj)\n",
    "\n",
    "# or could do it by hand: \n",
    "lj = (@analytic r -> r^(-12) - 2 * r^(-6)) * SplineCutoff(1.7, 2.5)\n",
    "set_calculator!(at, lj)\n",
    ";"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Users/ortner/.local/lib/python2.7/site-packages/matplotlib/backend_bases.py:2453: MatplotlibDeprecationWarning: Using default event loop until function specific to this GUI is implemented\n",
      "  warnings.warn(str, mplDeprecation)\n"
     ]
    }
   ],
   "source": [
    "# this should be done properly, using in particular \n",
    "# a proper scaling of the gaussian depending on temperature \n",
    "q = position_dofs(at) \n",
    "p = randn(length(q)) * 0.1\n",
    "set_momentum_dofs!(at, p) \n",
    "write(\"md.xyz\", at)\n",
    "# prepare for plotting\n",
    "pygui(true)\n",
    "PyPlot.ion()\n",
    "# create a counter to save only every 10 steps \n",
    "cnt = 0 \n",
    "dt = 0.003\n",
    "E0 = energy(at, q)\n",
    "Ep = [energy(at, q)-E0]\n",
    "Ek = [0.5 * norm(p)^2] \n",
    "for n = 1:1_000\n",
    "    p -= dt * gradient(at, q)    # Euler-A method\n",
    "    q += dt * p\n",
    "    cnt += 1\n",
    "    push!(Ep, energy(at, q) - E0)\n",
    "    push!(Ek, 0.5 * vecnorm(p)^2)\n",
    "    if cnt == 40 \n",
    "        cnt = 0 \n",
    "        # write current configuration into at\n",
    "        set_position_dofs!(at, q) \n",
    "        set_momentum_dofs!(at, p)\n",
    "        # write to file\n",
    "        write(\"md.xyz\", at, :append)\n",
    "        PyPlot.cla()\n",
    "        plot2d(at, \"n = $(n), Ep = $(Ep[end]), Ek = $(Ek[end])\") \n",
    "    end \n",
    "end \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjcAAAGgCAYAAABWo0bIAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAAPYQAAD2EBqD+naQAAIABJREFUeJzs3XdYlfX7wPH3Ye8hIuBgiIp7m6JmuVFz5CjNTM30W6lpZsNSo685sq/ltjJH7jJHZrm3uBA3KoqioAIOBGTDOc/vj/Pj5IkhIHAY9+u6nuvqnGec+xB47vP53M/9USmKoiCEEEIIUUYYGToAIYQQQojCJMmNEEIIIcoUSW6EEEIIUaZIciOEEEKIMkWSGyGEEEKUKZLcCCGEEKJMkeRGCCGEEGWKJDdCCCGEKFMkuRFCCCFEmSLJjRBCCCHKFEluhBBCCFGmmBg6gOKg0Wi4d+8etra2qFQqQ4cjhBBCiDxQFIUnT55QuXJljIzyPh5TLpKbe/fuUa1aNUOHIYQQQogCiIiIoGrVqnk+vlwkN7a2toD2h2NnZ2fgaIQQQgiRF/Hx8VSrVk33OZ5X5SK5yZyKsrOzk+RGCCGEKGXyW1IiBcVCCCGEKFMkuRFCCCFEmSLJjRBCCCHKlHJRc5MXarWa9PR0Q4dR5piammJsbGzoMIQQQpQjktwACQkJ3LlzB0VRDB1KmaNSqahatSo2NjaGDkUIIUQ5Ue6TG7VazZ07d7CyssLZ2Vma/BUiRVF48OABd+7coWbNmjKCI4QQoliU++QmPT0dRVFwdnbG0tLS0OGUOc7Ozty6dYv09HRJboQQQhQLKSj+fzJiUzTk5yqEEKK4SXIjhBBCiDJFkhshhBBClCmS3AghhBCiTJHkppQaNmwYKpUqy+bn51dsr9+nT59ieS0hhBAiP8r93VKlmZ+fHytWrNB7ztzc3EDRlD8hD0NYe3Etraq2onvN7oYORwghxP+T5OZfFAWSkgzz2lZWkJ+bi8zNzXF1dc3366hUKhYvXsy2bds4ePAgbm5uzJ49m/79++uOuXjxIuPGjeP48eNYWVnRr18/vvvuO2xsbPD39+eXX37RXQvgwIEDvPzyy/mOpbS6EXODF35+gfjUeACW91rO8CbDDRyVEEIIKOJpqZkzZ9KiRQtsbW2pVKkSffr0ISQkRO+YlJQURo8ejZOTEzY2NvTr14/o6Gi9Y8LDw+nRowdWVlZUqlSJjz/+mIyMjCKJOSkJbGwMsxVnUjVlyhT69evH+fPnGTx4MAMHDuTKlSsAJCYm0rVrVxwdHQkMDGTjxo3s3buXMWPGADBx4kRee+01/Pz8iIyMJDIyktatWxdf8CXApH2TdIkNwMQ9E4lNidU9vht/l6kHpvLzmZ9Ra9SGCFEIIcqtIk1uDh06xOjRozlx4gR79uwhPT2dLl26kJiYqDvmww8/5M8//2Tjxo0cOnSIe/fu0bdvX91+tVpNjx49SEtL49ixY/zyyy+sXLmSqVOnFmXopcL27duxsbHR22bMmJGncwcMGMA777xDrVq1mDZtGs2bN2fBggUArFu3jpSUFFatWkX9+vXp0KEDCxcuZPXq1URHR2NjY4OlpaVu5MjV1RUzM7OifKslyu3Y22y6sgmAM6POUNe5LjHJMcwOmA1AfGo8bVe0ZdrhaYz8cyTjd443ZLhCCFHuFOm01M6dO/Uer1y5kkqVKhEUFES7du2Ii4tj2bJlrFu3jg4dOgCwYsUK6tSpw4kTJ2jVqhW7d+/m8uXL7N27FxcXFxo3bsy0adP49NNP8ff3z/ZDNTU1ldTUVN3j+Pj4LMfkxMoKEhIK+Iafk5VV/o5v3749S5Ys0XuuQoUKeTrX19c3y+Nz584BcOXKFRo1aoS1tbVuf5s2bdBoNISEhODi4pK/QMuYRYGL0CgaOnh1oIlbE6Z3mM6rv77KglMLmNh6Ip/t/YxbsbcwVhmjVtQsDFzIiKYjaOza2NChCyFEuVCsd0vFxcUB/3wABwUFkZ6eTqdOnXTH1K5dG3d3d44fPw7A8ePHadCggd4HateuXYmPjyc4ODjb15k5cyb29va6rVq1anmOUaUCa2vDbPlt5mttbU2NGjX0trwmN6JgEtMSWXpmKQDjW2pHZHr79KaRSyMS0hJovay1bv++t/bxer3XAZhyYIphAhZCiHKo2JIbjUbD+PHjadOmDfXr1wcgKioKMzMzHBwc9I51cXEhKipKd8y/RwoyH2ce82+TJk0iLi5Ot0VERBT22yn1Tpw4keVxnTp1AKhTpw7nz5/Xmz4MCAjAyMgIHx8fAMzMzFCry18tyU9BPxGbEou3ozc9avUAtEXVMzvOBCDkkbambEyLMbzk+RL/bf9fALZf286NmBuGCVoIIcqZYktuRo8ezaVLl9iwYUORv5a5uTl2dnZ6W1mUmppKVFSU3vbw4cM8nbtx40aWL1/OtWvX+PLLLzl16pSuYHjw4MFYWFgwdOhQLl26xIEDBxg7dixDhgzRJZaenp5cuHCBkJAQHj58SHp6epG9z5IiOiGar498DcCnbT7FSPXPn0+3mt2Y02UOVe2qMrD+QGZ1mgVALada+NXQ9h764fQPxR+0EEKUQ8WS3IwZM4bt27dz4MABqlatqnve1dWVtLQ0YmNj9Y6Pjo7W3eLs6uqa5e6pzMcFuQ26LNm5cydubm56W9u2bfN07ldffcWGDRto2LAhq1atYv369dStWxcAKysrdu3aRUxMDC1atKB///507NiRhQsX6s4fOXIkPj4+NG/eHGdnZwICAorkPZYUl+5fotWyVsQkx1DPuV62t31P8J1AxIcRrO+3Hmuzf+qV3m32LgBrLq6RO6eEEKIYFGlBsaIojB07li1btnDw4EG8vLz09jdr1gxTU1P27dtHv379AAgJCSE8PFxX8Orr68v06dO5f/8+lSpVAmDPnj3Y2dnpPozLo5UrV7Jy5coCn1+5cmV2796d4/4GDRqwf//+HPc7Ozvnen5ZcuXBFdoub0tcahw1KtRg8+ubMTHK+59Ot5rdcLRwJCohisO3D9Peq30RRiuEEKJIR25Gjx7NmjVrWLduHba2trqpk+TkZADs7e0ZMWIEEyZM4MCBAwQFBTF8+HB8fX1p1aoVAF26dKFu3boMGTKE8+fPs2vXLiZPnszo0aOlG68oFmN2jCEuNY5WVVtxYsQJajnVytf5ZsZm9KujTd5/Df61KEIUQgjxlCJNbpYsWUJcXBwvv/yy3tTJr7/+8w/8999/zyuvvEK/fv1o164drq6ubN68Wbff2NiY7du3Y2xsjK+vL2+++SZvvfUW//3vf4sy9FJr7dq1WXrfZG716tUzdHhFbvu17TT9sSmT9k5CUZTnvl5oTCj7w/ZjrDJmQ78NOFk5Feg6fetoezf9ff3vQolLCCFEzop8WupZLCwsWLRoEYsWLcrxGA8PD/7+++/CDK3M6tWrFy1btsx2n6mpKZC3/y+l0b0n9xj4+0AS0xM5G3WWJm5NeK3ea891zT+u/gHAy54v4+HgUeDrvOT5EubG5kTER3Dl4RXqOpffKVVRNt2Jv8P8k/MxMTJhgu8EKlpVNHRIohyTtaXKGFtbW2xtbQ0dhkEsDVpKYvo/t6//GPTjcyc3W0O2AtCn9vOtgG5lasXLni+z68YudlzfkWNyoygKoTGhXHt0jdiUWBq5NqJ+pfrP9dpCFLXIJ5G8sPQFIhMiAdgWso3AkYFYmloaODJRXhVrEz8hitKf1/4EYPKLkwE4fPuw3vpP+ZWcnsyJO9p+QK/UeuW54+vq3RWAvWF7s92/P2w/tRfVptbCWryy/hXe3PImDZY04P2/3i+zo22i9FMUhaFbhxKZEEl1x+pUtKpI8INgfgz60dChiXJMkhtRJsSlxHE26iwA7zZ/Fw97DzI0GQTeDSzwNS9EXyBDk4GzlTMe9gWfksrUqbq2E/fh24dJU6fp7dt3cx9dVnfh2qNrmBub09ClIW3d26JCxZLTS1gcuPi5X1+IonAk/Ah7bu7BzNiMv9/4m6/ba3tBzTs5jwxN0SxwLMSzSHIjyoSAiAA0igZvR2+q2FWhZVVt3dGpu6cKfM2gyCAAmlVuhiq/a2Nko36l+lSyrkRSepJuRAggIS2Bt7e9jVpRM6DuAO5/fJ/z757nyPAjfN/1ewCmHZ5Gcnryc8cgRGGbeVTbnXt44+H4VPRhSKMhOFk6cSv2FrtCdxk4OlFeSXIjyoRDtw4B8JLHSwA0dW0KwKUHlwp8zaB72uSmuVvz54xOS6VS6UZv9t78Z2rqp6CfCI8Lx9PBk+W9l2Nn/k9H7fdavIe7vTvRidH8cv6XQolDiMJyJvIMO0N3YqQy4pM2nwDa+rI3GrwBwLpL6wwZnijHJLkRZcKh2/+f3Hhqkxufito1sEIehhT4mk+P3BSWjl4dgX+SmzR1Gt+f0I7OfPHiF9iY2egdb2ZsxriW4wBYc2FNocUhRGGYdVS7zMjA+gOp7lhd9/zgBoMB2Hp1KwlpCQaJTZRvktyUQyqViq1bt+a439PTk7lz5xZjRM8nIS2B0/dOA/+M3NSuWBvQLmRZkGLc5PRkLt3Xjvo0cyu85CZz5ObU3VPEpcSx4dIG7sTfwdXGlTcbvpntOa/Xex0VKgIiAoiIk0VgRdFSa9TMDpjNW1ve4kzkmRyPuxh9kd8v/w7AZ20+09v3QpUX8Hb0Jik9SddOQYjiJMlNKTVs2DBUKlWWzc/P77mvHRgYyKhRowohyuJxLOIYakWNh72HrhdNdcfqGKuMSUhL4N6Te1nOSUxL5FbsrRyveSH6AmpFjbOVM1XtquZ4XH6527tTs0JN1IqanaE7mX5kOgDjWo7DwsQi23Oq2FWhrbt2zbDMDxMhisp/D/2XT/d+yuoLq+m0qlO2fz8ZmgzG7BiDgkL/uv1p4NJAb79KpSrxU1NHbh9h6NahfBvwrRQ+l0GS3JRifn5+REZG6m3r169/7us6OztjZWVVCBEWD129zf9PSYF2Ose7gjegHb152vmo83jM9cBrnhdvbHoj28UsC7uY+GmZ3YoHbhrItUfXcLJ0YnSL0bmek7l8wx8h8i1YFJ2bj2/yTcA3usePUx7z9eGvdY8PhB1g6Nah1FxQk8O3D2NpYsmcLnOyvVZmcrMrdBcPEh8UbeD5FBAeQIdVHVh1fhWf7P2E8TvHGzokUcgkufkXRVFITEs0yJbf6RNzc3NcXV31NkdHx3y/5y+//BI3NzcuXLgAZJ2WUqlU/Pzzz7z66qtYWVlRs2ZNtm3blu/XKSqZ9Tbt3NvpPe/jpK27ufrwqu45RVEY+edIHiU/AmD9pfW6mpenFXYx8dPeb/E+Vqb/JI+T203G1jz3xou9a/cGtLfdPkp6VOgxCaEoCh/u+pBUdSqdqnfi0DDt39XPZ37mduxtlp1ZpksIbsXewtbMll/7/4q7vXu216tdsTZN3ZqiVtRsvLyxON9KrpLSkxj2xzAyNBnUc9YuSbMocBFHw48aODJRmKRD8b8kpSdhM9Pm2QcWgYRJCVibWRfb6ymKwgcffMD27ds5cuQINWrUyPHYr776itmzZ/Ptt9+yYMECBg8ezO3bt6lQoUKxxZudpPQk3e3eT4/cALoFLq8/uq57LigyiMB7gViYWDCl3RS+2P8FX+z/Ar8afnqdgE/c1d6q3bxy4Sc37vbubH5tMwsDF9KySktdwXBuPB08aejSkAvRF/jr+l+81eitQo9LlG8rzq1gW8g2TIxMmNt1LvUq1aNT9U7svbmXF35+QTf68kaDN3izwZu0rtYaewv7XK85uMFgzkSeYfWF1bzf4v3ieBvP9OmeTwmNCaWKbRWOvn2Uj3Z9xPJzy5m8fzIHhx00dHiikMjITSm2ffv2LItjzpgxI0/nZmRk8Oabb7Jv3z6OHj2aa2ID2hqfQYMGUaNGDWbMmEFCQgKnTuW/h0xUQhQjt43Ea54XLX9uyZLAJaRkpOT7OplO3DlBuiadyraV8Xb01ttXs0JNAK7H/JPc7AzdCUD3mt2Z1HYSPWv1JE2dxltb3iJdnQ7Ao6RHXH5wGYA27m0KHFtuutboyp+D/mRyu8l5nvbq7aMdvZGpKVHYrj26xgc7PgBgWvtp1KukHdGY0WEGxipj7ifeR0Hh/ebvs+bVNXSr2e2ZiQ3AoPqDMDM248SdE3rtDwwh8kkkswNmszBwIQDLei3DwcKBr9p/hamRKYduH5LRmzJERm7+xcrUioRJhrl18empirxo3749S5Ys0XsuryMpH374Iebm5pw4cYKKFZ+9wF3Dhg11/21tbY2dnR3379/PV7xRCVG8uOJFQmNCAbgVe4tTd08x/ch0vnjxC95u8jbmJub5uua+m/sA7cKW/04SajplTW5239gNQJfqXVCpVPz4yo8cXXyUs1Fn+WL/F8zuPJuAiAAA6lSsU6IW/+tTuw/TDk9jV+guUjJScixAFiI/0tRpvLHpDRLTE2nv2Z6PW3+s29eiSgu2DtzK0jNLaVutLRNbT8xXDZqbrRvvNX+PeSfn8dnezzg18hRGquL9Tp2QlsCQLUPYevWfO0Q/8v2IrjW0y6FUtavKsMbDWHpmKdOPTGfH4B3FGp8oGjJy8y8qlQprM2uDbPktXLW2tqZGjRp6W16Tm86dO3P37l127cpbB9HMFcWf/jlpNJo8x6pRNAzePJjQmFA87D3YNnAbc7vOpapdVe4+ucv7f79PrYW1WBq0VDeCkhe7bmjjz1y36WmZIzc3H98kQ5NBfGo8x+8cB6CLdxdA+4/v0p5LAfj22LdsC9nG4duHAXjR/cU8x1Ecmrg2oZpdNRLTE3VJnRDPa9bRWQRFBuFo4ciqV1dhbGSst/+VWq/wx8A/+LjNxwUqrs/s3xQUGWSQu/0+2PGBLrFpXrk5c7rMYXbn2XrHfNb2M4xVxuwM3alrKyFKN0luyqlevXqxbt063nnnHTZs2FDkr7fu4jr2h+3HytSKXW/uoqdPT8a1Gkfo2FAWdluIm40b4XHhjNo+Cp+FPvxw+odnLnr5IPGBrg9HZrLytCp2VbAwsSBDk8Ht2NscCDtAhiaDmhVq4uXopTuuX91+urqX1za+xpzj2rs/Ont3Lqy3XyhUKhW9fHoB8GvwrwaORpQE8anxTN4/mf6/9Wf9xfzfKfk4+bHu931R90WF2vYgk7O1MxN9JwLwxf4v8vXl5XmFPQ5j5bmVAOx/az+BIwOZ4Dshy+hRdcfquru7ZhzJ29S+KNkkuSnFUlNTiYqK0tsePnyY5/NfffVVVq9ezfDhw/n996L7RpWhyeC7498BMLvTbF33YABzE3NGvzCaGx/cYG7XubhYuxAWG8Z7f71H1e+qsunyphyvu+fmHhQUGrk0wtXGNct+I5WRrg4nNCZUN8qTXSI0u/NsWlVtRao6FQAXaxd61upZ8DddRIY0HAJok8WnC6VF+ZOUnkSnVZ2YfmQ6m65s4o3Nb7Dg5IJ8XWPeyXnEp8ZTv1J9Xq//ehFFChN8J+Bs5UxoTCibr2wustf5t/WX1qOg0NGrI+292ud67KS2k1ChYsvVLboGnqL0kuSmFNu5cydubm56W9u2bfN1jf79+/PLL78wZMgQNm8umn90YlNiSUpPwreqL++1eC/bYyxNLRnXahw3x93k+67f4+Pkw5O0JwzYOCDHKZjckpVMmYlU8INgXb1NdlNYZsZm/D7gd7p4d6Guc11+G/Bbvut/ikPLqi3pXrM7akXN67+/zpYrW545wiXKpo93f0zgvUCcLJ14u/HbAHyy95M8LzkSmxLL3BPalg9T200t0loYW3NbXS+neSfnFdnr/Nv6S9rRrEH1Bz3z2DrOdehXV9tP6uleP6J0UikF6U1fysTHx2Nvb09cXBx2dnZ6+1JSUggLC8PLywsLCynQLGyP4h8RGBzIewHvsWnwJpq6Nc3TeRmaDN7+421WX1iNh70Hl0df1iu4VmvUuM5x5WHSQ/a9tY8OXh2yvc60Q9OYenAqrau15ljEMUyMTIj5JOaZfWVKsmuPruG7zJeY5BgAXSO1nBJHUfbciLlB7UW1ydBksGfIHjp6dcRvrR+7b+ymi3cXdg7e+cz6mMy/jbrOdbn43sUiL/SNSojC/Xt30jXpnHrnFC2qtCjS17sVewuveV4Yq4x58PEDHC2f3QPs9L3TtFjaAlMjUyI+jMDFxqVIYxTPltvnd25k5EYUGUVRiEyIBGBAvQF5TmwATIxMWNxjMe727tyOu82Ksyv09h+LOMbDpIc4WDjkWvib+ZrHIo4B0Lpa61Kd2IC2f0/QqCDGtRyHt6M3yRnJvP/3+/wZ8qehQxPFZO6JuWRoMvCr4Uen6p1QqVQs7r4YM2Mzdt/YzdqLa3M9Pz41Xte8ckq7KcVyB5Orjatu6iu7xpmFLbNzeYsqLfKU2IC24LhV1Vaka9JZemZpUYYnipgkN2XM2rVrs/S+ydzq1atXrLHEJMeQkpGCkcqI8S3z397cxsyGT9t8CmjvZHq6EHFbiLZDco+aPTA1Ns32fIAXPV7ExOifjgd9fPrkO46SyNPBk7l+c7k+9jpjXxgLwMg/RxKXEmfgyERRU2v+6fj7wQsf6J73ruCtW8ByyJYhOMxyoPWy1pyPOp/lGjOPzORxymNqV6zNgLoDiidw4MNWHwLwW/BvhMeFA9o7KYuiyDizc3nmYrp5NabFGACWnF4ia06VYpLclDG9evXi3Llz2W5///13scWhKIpuwT07czsqWBWsk/HwxsNxtnLmdtxtNlzaoLt2ZiO7zMZ2ObEzt6N/3f4AOFo45rjydmmlUqmY3Xk2Pk4+RCdG8+2xb/N8rkbRkJiWWITRlR7fHf8Ox28cabCkAReiLxg6nFwdvn2Y6MRoHC0c6Vi9o96+KS9N4Z0m76BCRVxqHMfvHKfDqg66ppSgXRT2uxPaAv9vOn2T5dbvotTUrSkdvDqgVtSM/HMkE3dPpOLsitjPsmfGkRn5XoImNwdvHQS0PbDyo3/d/jhbOXPvyT3+uvZXocUjipckN2WMra1tlt43mZuHh0exxfEo+RGp6lSMVcbPNQ1kaWqp+7Y34+gMNIqGM5FnuB5zHTNjM/xqPHsV9J97/sx8v/kcffsoztbOBY6lpLIwsWBWp1mA9kM68knkM8/JXDzUdqYtH+z4oFA/VEqbA2EH+Gj3R8SmxHLp/iX6bOhDQpphGnnmRWYbgL51+mJmbKa3z8TIhKW9lvLok0ecf/c8LSq3ICY5hj4b+hCbEkt4XDh9NvQhTZ1Gj5o9DHJH4OxOszExMmH3jd3MOT6HxymPSc5I5ov9X/DJnk8K5XcxIi6CsNgwjFXGtKmWvy7j5ibmDGs8DICfzvz03LEIw5DkRhQ6RVGITogGwMnK6bnn80e/MBoHCweuPrzKpsub+DHoR0D7DSsviZO1mTVjW46lrnPd54qjJOvt0xvfqr4kZyTjf9A/12PT1em8/vvr3Im/g4LCglMLdHfNlDeKovDZPu1UTm+f3lSzq0ZYbNgzf4aGkqHJYNMVbXuE1+vlfOu2o6UjDV0a8tcbf+Fu7871mOu0Xd6W5j81Jyw2jOqO1Vn16qpCX/E+L5pVbsa2gdto59GO3j69+WPgH3zfVVuD87/j/+Oj3R+Rpk4DtO83KiGK4PvBHLl9hJuPb+bpNTKnpJq6NS3Ql6t3mr4DaJdryZw+Kw/UGnWZuftSkhtR6BLSEkjOSMZIZUQFy+dfWNPO3E7XZG/knyNZdnYZAKOajnrua5cVmdNTAMvOLtNbCf3fNl7eSMijEJytnPF/yR+Az/Z9xpUHV4oj1BLlXNQ5Tt09hYWJBT++8iM/vPIDAAtPLSzUD7XIJ5Hsu7nvmTVRt2JvseLsCvaH7UetUWfZfyDsAA+THlLRquIz+7aAtoHe1te3YmtmS/CDYB4kPaCecz0ODj1YKH+bBdWtZjcODTvE1oFb6eXTi/GtxrOo+yJAW2xc4ZsKVPimAqbTTHGb40b9JfVpt7Id3vO96bK6C4+SHuV6/YJOSWWq5VSLlz1fRqNoWH52ebbHKIrCnht72HBpA6kZqQV6nZIk8kkkdRbVoeLsijm+59JEkhtR6O4natecqmBZQa+Y93mMazmOyraViUuNQ6NoeNH9Rdp5tCuUa5cVbd3b0sunl66eIadvYIsDFwMw9oWxTH1pKj1q9iBNncao7aPQKHlfUqMsyKzd8qvhh4uNC91qdONlz5dJVafy1cGvnvv6GkXD7IDZeMz1oNPqTnjP986xvf/O0J3UXlibt7e9TcdVHWn8Y+Ms9T+ZU1L96vTL899WE7cmXHjvAjM7zmRF7xWc+c8ZqtlXe743VgTeb6FdlNPVxpXE9EQepzwGQIUKJ0snvBy0t3XvubmH3ht651rsW9Bi4qdlfnladnZZtonm1ANT6bKmC4M2DaL7uu6lvvj4y4Nfcj3mOumadMbuGMvDpLw3hC2JJLkRhSpDk0FsSiwAlawrFdp1HS0d2fXmLnrU7MHA+gP5bcBvBhlSL+m+6fQN1qbWHA0/SoMlDXR3lWU6H3WegIgATIxMeKfpO6hUKhZ1X6Q75+czPxsocsP4d2G6SqViZseZAKw8v7JAxcVqjZqDtw4y/fB0Wixtwad7PyVdk46duR2Pkh8xaNMgktKT9M65FXuL1za+Rqo6lbrOdbE3t+fS/Uu0XtaaA2EHAHiS+oTfgn8DYGD9gfmKydPBk8/afsawxsOy1OmUJIMbDub2+Ntcfv8ywe8Hc3/ifdKnpPPwk4fcHHeTc++ew87cjoCIgByXm7gbf5fQmFCMVEa0dc9fU9OnvVrnVZwsnbgTf4edoTv19gXfD2b6kem6x/vD9rPq/KoCv5ahZWgy9Nb9SkpPyrU7fGkgyY0oVI+TH6OgYGlime9Vzp+lfqX6bH9jO+v7rc92uQUBtSvW5uCwg3jYexAeF07vDb0Z+PtAktOTAVhwStuev2+dvrjZugHg4eDB1x2+BuCTPZ/kqSC5LLgde5tzUecwUhnxSq1XdM+3qtqKvnX6olE09P21LxeiL+S5yDU2JRbfZb60/6U9kw9M5kzkGWzNbPm558+Ejw+nql1VQmNC+eaofgdc/4P+PEl7QutqrTlI4mF9AAAgAElEQVT7n7Pc+OAGHb06kpieyCvrX+HQrUOsvbiWJ2lPqOVU67lGJEo6M2Mz6jjXoa5zXZytnfXu5qpfqb7udvfpR6ZnO6KSOWrTxLUJ9hb2BY7DwsSCtxq9BcDi04v19n177FsUFPrW6cv/Ov8P0C5Aml08pcHxiOM8TnlMBcsKuuR+Q3DRrzlYlCS5EYXqUbJ2LtzJysnAkZRfzSs3J/j9YD5r8xkmRib8GvwrXdd05fKDy7rmbpm9cTKNfWEszSs3Jy41jgEbB+gKwsuyzFGttu5tqWhVUW/fDz1+wN3enRuPb9Doh0bYzLShzfI2uhXjczJ2x1gC7wViY2bDwPoD+bbzt1wfe50RTUdgb2GvK5z9/sT3ug7TVx9eZfWF1QDM7ToXM2MznKyc2P7Gdvxq+JGUnoTfWj8m7tYuPvlus3fL9ajl6BdGY29uT8ijEI6GH82yf3/YfgDaez67JulZ3mv+HipU/H39b916Uw+THuraUnzc+mP+0/w/VLCswPWY6+wI3fHcr2kIf1/Xtgnxq+GnW6ri0K1DuhKD0kiSm1Jq2LBhqFSqLJuf37NvjVapVGzdujXfr/ms89Iy0nS30FawMFyxotDeITaz00z2DtmLnbkdR8KPUG9xPVIyUvCt6pvl9lhjI2OW91qOvbk9AREBVP2+Kj4LfZiwa0KWKZSyYmuI9nc5u15JztbOHBp2iF4+vTBWGZOUnsSxiGN0Xt2Zk3dOZnu9wLuBrLmwBiOVEbvf3M36fuuZ2HqiXgv/vnX60silEU/SnjDnmHY17qkHpqJRNPT26a23JIGFiQVbXt9CV++upGSkkJieSD3nerzb/N3C/DGUOnbmdvSt0xdAN033tAO3tNN4eSm4fpaaTjV1603NDtAW7P8U9BOp6lSaujWlZZWW2JjZMLTRUABWnFuR47VKsmN3tB3cO3l1wsPBgyauTVBQ2BW6y8CRFZwkN6WYn58fkZGRetv69dnPQxeH2FRtrY2NmQ1mJiV3Xr88ecnzJQ4PO6ybxnOydOLnXj9n+82/gUsD9g/dTzO3ZmRoMrj26Brfn/iezqs756vvi0bRsDN0J6vPr9bVX5U0j5Mf69rz59QI0tPBkz8G/kHi54mEjAmhZ62epKnT+M/2/2RbPPrDae2dVgPrD8S3mm+21zRSGfHVy9pC5fmn5vNnyJ9svLwRFSqmtZ+W5XgLEwv+GPgHP77yI7M7zebw8MNYmloW6D2XJa/Vew2A36/8rjcVdDv2Njcf38RYZZzrsiz5kdklfd3FdZy8c1K38Of4luN1f0fDGw8H4M+QP3mQ+KBQXre4aBQNZyPPAtrb9AG61+wOwN+hxdf4tbAVaXJz+PBhevbsSeXKlbP91p/d6MO/Rx5iYmIYPHgwdnZ2ODg4MGLECBISirDBlqJAYqJhtnw2rzI3N8fV1VVvc3TMfQ0VT09PAF599VVUKpXuMcCSJUvw9vbGzMwMHx8fVq9enafzMmV+kDlYOOTrfYii1ci1EZffv8yOwTsIGROSa7+fpm5NCRwZSPj4cH7t/ysOFg4cizjGkC1D8nQnVWpGKr039Kbb2m68tfUt6i+uz/VH1wvz7RSKv6//jVpRU8+5Ht4VvHM91tzEnFpOtVjeezkVLCtwPvo8ay6s0TsmNiVWV6PwXvPcFzDt5dOLZm7NSEhLoNeGXgAMbTyUBi4Ncnz9Uc1G8XGbjw16+3ZJ0tGrI3bmdtxPvK9391nmqE3zys0LbQ255pWb069OP9SKmlbLWnE/8T7V7KrpFXU3cGlAM7dmpGvSs13XS1EUVp1fRb/f+jHr6KwSdWfVzcc3eZL2BHNjc+pUrANAtxrdANgVuqtExZofRZrcJCYm0qhRIxYtWpTjMf8effj3yMPgwYMJDg5mz549bN++ncOHDzNqVBH2N0lKAhsbw2xJRT/8HxgYCMCKFSuIjIzUPd6yZQvjxo3jo48+4tKlS/znP/9h+PDhHDhwINfzMqk1ap6kPgHA3rzgRXyiaDhaOuJXwy9PtVAqlYpq9tV4rd5r/PXGX5gZm7H16lam7J+S63lqjZpBmwax/dp2LEwsqGpXlbtP7vLmljef+Q/k7hu7eXf7u3x9+GtdLUpRym1KKicVrSrqvsV/E/CNXrK35sIaktKTqOdc75kdcVUqFSv7rNTV+dR1rqsrShV5Y2psShfvLsA/9SLw1JRUIdTbPG1h94W4WP8zvTjPb16WNe3ebvI2AMvPLs9SgD71wFSGbh3K5iubmbRvEh/s+ICSInPUpoFLA917alm1JY4WjjxOeZzjNGyJpxQTQNmyZYvec0OHDlV69+6d4zmXL19WACUwMFD33I4dOxSVSqXcvXs3z68dFxenAEpcXFyWfcnJycrly5eV5ORk7RMJCYqiHUMp/i0hIc/vaejQoYqxsbFibW2tt02fPv2Z52b3/6J169bKyJEj9Z4bMGCA0r1791zPyxSTFKME3g1ULkRdUDQaje75LD9fUeqsOrdKwR8Ff5TXN76uzD0+V1l5dqWy7eo25W689u8wQ52hvLf9PQV/FLNpZsreG3uVO3F3FPuZ9gr+KItOLcrx+vNPzNddH3+U6vOqK6GPQovs/TxMfKiYTzNX8EcJuheUr3PjUuIUh1kOCv4omy5vUhRFUTQajVJvUT0Ff5QFJxfk+VqPkx8rR28fVZLT5W+jIJafWa7gj9L8p+aKomj/P1T7rpqCP8ru0N2F/np34u4oc4/PVY6FH8t2f0xSjO736vTd07rnvz/+ve53e9Dvg3T/HXw/uNBjLIhJeycp+KOM3Kb/7//A3wcq+KN8vvdzA0Wmldvnd24MXnNz8OBBKlWqhI+PD++99x6PHv3TefL48eM4ODjQvHlz3XOdOnXCyMiIkydzziZTU1OJj4/X2/LMygoSEgyzWeXv1un27dtnWRzz3XcLVmx45coV2rTR/8bZpk0brlzJW9fauFRt51V7C/tyfSdHWTSk0RC+fOlLQNtEbvyu8Qz7Yxi9NvSiyndV8JrnRa2FtVhyegkAq19dTcfqHaliV0V3i/k3Ad/oWuo/LeheEBN2T9C+TsMheDp4cvPxTXpt6FVki3quOr+KVHUqTVyb0MS1Sb7OtTO3060anbnQY0BEAMEPgrEytWJIwyF5vpaDhQNt3NtgYWKRrxiEVua6cqfvnSY6IZprj64RER+BqZEprau1LvTXq2JXhXGtxuVYT+Vo6cirdV4F/iksXhq0lA93adfGm95hOuv6raNfHW2B8jcB32R7neJ2Nko7cvPvv4XuNbR1N1tDtpbKtecKp31sAfn5+dG3b1+8vLy4ceMGn3/+Od26deP48eMYGxsTFRVFpUr6jeBMTEyoUKECUVFROV535syZfPVVAbuLqlRgbV2wc4uZtbU1NWrUMHQYADIlVcb5v+yPXw0/fr/8OxHxEcSlxHHvyT2CHwRzK/YWoC0k/6HHD7piT9Cu0TP9yHTC48JZc2GNbugetHUIH+76kAxNBv3r9ueXPr8QlRBFs5+acfnBZT7c9SE/9Xy+hQtjU2I5cvsIsSmxuNu7Y2ZsxqwA7SKjo5qNKlAiPq7VOL478R1BkUHsubmH745rV9ge3GDwc/VVEfnjZutGU7emnIk8w64bu4hK0H4mtPdqj7WZYf4Nf7vx22y4tIHlZ5eTrk5n6ZmlAEz0nciktpMAbYHypiubWH9xPXO6zMnShuB5KYrCktNL+DX4V3ycfJjWfpreHXv/ljkt1dStqd7zr9R6BStTKy4/uMyBWwfo4NVBr5Hka/VeK7S6pqJg0ORm4MCnCrIaNKBhw4Z4e3tz8OBBOnbsWODrTpo0iQkTJugex8fHU61ayWs3biimpqao1frNpurUqUNAQABDhw7VPRcQEEDdunVzPQ+0RaSpau3aKjZmNkUUtTC0VlVb0apqK73n4lPjOXnnJAlpCbTzaJelpsfCxIKJvhOZuGcis47OYmijobqmbHtv7uVI+BHMjc2Z23UuKpUKN1s31vVbR/tf2rP0zFL61umbp5Xf/y0qIYpJ+yax/uJ63e/m0+o519NLtPKjolVFRjYdybyT8+j7a18S0xMxMTLhs7afFeh6ouC61+jOmcgz/HX9L+7E3wHyV0dV2DpV70Tn6p3Zc3OPbkXxsS+MZXbn2bpEukWVFjRza0ZQZBBrL6xlXKtxhRrDrKOz+Hz/5wAcvn2Yw7cPEzQqKNuELyohiujEaIxURlkK2h0tHRnWaBiLTy/G/6A/iWmJfLDzA92Xme9OfMfuN3dTxa4KZyLPoEJFE7f8jYQWJYNPSz2tevXqVKxYkdDQUABcXV25f1+/iVBGRgYxMTG4uubcodbc3Bw7Ozu9rSxKTU0lKipKb3v48NnrgXh6erJv3z6ioqJ4/Fi7fsvHH3/MypUrWbJkCdevX+e7775j8+bNTJw4MdfzAJ6kaUdtrE2t9bqJirLPztyOzt6dta3qcyhWfrrJWeaK1oqiMOWAtkD5vebvUcWuiu74lz1f1i2UOmLbCB4nP8560VyEPAyh6Y9NWXluJanqVHycfOjg1QFvR2/szO3oVqMbOwbveK5lCD7y/QhLE0sS07VTZ6NbjKa6Y/UCX08UTGZn6d+Cf+NYhLZXSy+fXgaLR6VSsem1TXzk+xG9fHrx+4Dfmec3L8sIoa74+FzW4uPncSv2Fl8d0s5afOT7EVXtqhLyKISpB6Zme/z5qPOAdqHQ7DrKT2w9EWtTa46EH6HXhl7cir2Fu707zlbOXH5wmU6rOxGdEM3wP4bTfGlzfr30a6G9l+dW6NU/OSCXYtRMERERikqlUv744w9FUf4pKD59+p/irF27dhVtQXEpMXToUAXIsvn4+Dzz3G3btik1atRQTExMFA8PD93zixcvVqpXr66YmpoqtWrVUlatWpWn827G3FQC7wYqd+LuZHmt0vrzFYXrywNfKvijNP6hsaLWqJXNlzcr+KNYTbdSop5EZTk+MS1RqbWgloI/yuBNgxWNRqOcjzqvzDk2R9lyZYui1qizfZ2ktCSl5vyaCv4odRfVVY6FH9MrcC9Mx8KPKYN+H6RMPzxdSctIK5LXELnTaDRK4x8a64p0e67raeiQ8uTp4uOcitrT1enK9UfX8/W7NXzrcAV/lA6/dFA0Go3y97W/FfxRTP5roiv+f9qsI7N0NwrkZOf1nUqVOVUUs2lmyui/RiuxybFK2OMwpep3VRX8UYy/MlbwR6nwTQXlQeKDPMeaVwUtKC7S5ObJkyfK2bNnlbNnzyqA8t133ylnz55Vbt++rTx58kSZOHGicvz4cSUsLEzZu3ev0rRpU6VmzZpKSkqK7hp+fn5KkyZNlJMnTypHjx5VatasqQwaNChfcZTF5KakyPzQCbwbqMSlyM9XZO9h4kPFerq1gj/KpL2TdHe1fLHvixzPOR5xXDH6ykjBH8XlWxe9O6pe3/h6tgnOlP1TFPxRKs+prNxPuF+Ub0mUEGcjzyo+C3wUnwU+StjjMEOHk2eZdyON/mt0ln2Xoi8p3vO8FfxRvOZ6KZfvX37m9R4lPVIsvrZQ8EcJCA/QPd9mWRsFf5Sp+6dmOSfz7q0Zh2fkem2NRqNkqDP0nrt8/7LiPNtZlzxtuZL74EVBlci7pU6fPk2TJk1o0kQ7DzdhwgSaNGnC1KlTMTY25sKFC/Tq1YtatWoxYsQImjVrxpEjRzA3N9ddY+3atdSuXZuOHTvSvXt32rZty08/PV+RoSg86ep03V0w1qaloxBbFD8nKydmd9a2r595dCYR8RF4Onjy+Yuf53hOq6qtmO83HxMjE6ITozExMqGLdxdMjUz5NfhXVpzVb3Ufnxqv6x47t+tcnK2di+4NiRKjsWtjro65ytUxV/F08DR0OHn2dmPt1NTai2tJyUjRPf8w6SE91vXgxuMbAITFhtHvt36kZmStHXvaL+d+ISUjhUYujfCt+s8dXZnryC09szTLwp7no7XTUo1cG+V6bZVKlaXkoI5zHa6MvsKGfhsIfj+YPrX75HqNYlckqVYJU55GbtasWZOl903mVrdu3UJ/vUdJj5TAu4E59mwoaz9fUXAajUb54O8PFOOvjJUa82so5yLP5em8O3F3lL039iqPkx8riqIo/wv4n4I/itM3TrrnFEVR5hybo+CPUnth7RynrYQoKTLUGYr79+4K/ijrLqxTFEU7FdXhlw4K/ig15tdQQh6GKJW+raTgjzL9cM49zDQajW4a94fAH/T2pWakKo6zHBX8Ufbd3Kd7Pjk9WTellF1JQUlRIkduRPHr1atXlt43mdvffxf+OiGZvUhk1EY8i0qlYl63eSR9kcS1Mdee+W0xUxW7KnSs3lG3rMcHLT+grnNdHiU/4tuAbwHtCGLmqM1Hvh9hpJJ/2kTJZmxkrFuTyv+QP8npyYzfOZ79YfuxMbNh6+tbqeVUizldtAuszjo6i4dJ2d8wsj9sP9ceXcPGzIY3Gryht8/M2Iz+dfsD2vWxMgXfD0atqHGydKKybeWieIsGJf8ClDG2trbUqFEj283Dw6PQXy9zQUW5BVzklZmx2XM1ejQ1NmVGhxkAzD05l8gnkay9uJbwuHAqWVfizYZvFlaoQhSp8a3G42bjxrVH16g+vzqLArVLFf3S5xfqVaoHwBsN3qCxa2OepD1h1tFZ2V5n7sm5ALzV8K1se89kJjy/X/5dN72VuSZXY9fGZbLxqiQ3/08phR0YDU2jaEhK166HlVPTLPm5iqLQy6cXvlV9SUpPYtzOcXx9WNsJeaLvROn4K0oNBwsHVvRegbmxOVEJURipjPjplZ/oW6ev7hgjlZEumV94aqGun0+myw8us/3adlSoGN9qfLav086jHVVsqxCXGseO0B0AHLp9CKDQVk8vacp9cmNsrC2SSkvL2hpe5C4pPQkFBRMjE8yNzbM9JvPnmvlzFqIwqFQqZnXSfovdeHkjNx7fwM3Gjfda5L4itxAlTdcaXQl+P5gVvVdw6b1LjGw2MssxfjX8eNH9RVLVqXx1UL/7/sTd2l5kfWr3oaZTzWxfw0hlxKD6gwDt1JSiKBy8dRDQ9pUqi1RKOfhqHR8fj729PXFxcVka+imKQnh4OOnp6VSuXBkjo3Kf7+XZw8SHRCVGYWNqg6ejZ5b9Go2Ge/fuYWpqiru7e5kc+hSGNffEXCbtm0QFywr8PuD3HNf9EaK0CwgPoO2KthirjDnznzM0dGnI8rPLGbFtBKZGplx6/xK1nGrleP7ZyLM0/akpFiYWHBl+hBZLW2BubE7sZ7ElerQzt8/v3JT75Aa0owthYWFoNBoDRFd6PUh6QFJaEg4WDjmuqWNkZISXlxdmZgXvBitEbtLV6RgbGUsRsSjz+mzowx8hf+Dp4Mnr9V5nzvE5ZGgymNZ+GpPbTc71XEVRqLu4LlcfXqWaXTUi4iPoXL0zu4fsLqboC0aSm1zk5Yej0Whkaiqf/Nb4cSv2Fst6LaONe5tsjzEzM5PRMCGEKASPkh7R8ueWuh44AG81eosVvVfkKbmfeWSmbt0pgPX91jOw/sBczjC8giY3Bl04syQxMjLCwqLkDs2VNAlpCRy+exgFhYZVGsrPTgghipiTlRMBbwcw48gMbsXdom/tvrzV6K08T/mPazWObde2ceLOCXrW6slr9V4r4ogNR5IbUSAXoy+ioFDZtrJ0ghVCiGLiYuPCvG7zCnSulakVR4cf5VbsLao7Vi/TdZCS3IgCORd1DtD2SBBCCFE6GBsZ413B29BhFDkphhAFkpncNHLJW5dZIYQQorhIciMKJHPBNRm5EUIIUdJIciPyTa1RcyH6AiDJjRBCiJJHkhuRb9djrpOckYyVqRXejmV/7lYIIUTpIsmNyLfzUdopqYYuDTE2kmUVhBBClCyS3Ih8090p5SJTUkIIIUoeSW5Evp2L/v87pVzlTikhhBAljyQ3It8yp6WkmFgIIURJJMmNyJfohGgiEyJRoaJBpQaGDkcIIYTIQpIbkS+Z/W1qOtXE2szawNEIIYQQWUlyI/JFpqSEEEKUdJLciHzRFRPLsgtCCCFKKEluRL7IgplCCCFKOkluRJ4lpycT8jAEkORGCCFEySXJjciz4AfBqBU1Fa0q4mbjZuhwhBBCiGxJciPy7OkpKZVKZeBohBBCiOxJciPyTHenlCy7IIQQogST5EbkmSy7IIQQojSQ5EbkiUbRSI8bIYQQpYIkNyJPbsXe4knaE8yNzfFx8jF0OEIIIUSOijS5OXz4MD179qRy5cqoVCq2bt2qt19RFKZOnYqbmxuWlpZ06tSJ69ev6x0TExPD4MGDsbOzw8HBgREjRpCQkFCUYYtsZBYT16tUD1NjUwNHI4QQQuSsSJObxMREGjVqxKJFi7LdP3v2bObPn88PP/zAyZMnsba2pmvXrqSkpOiOGTx4MMHBwezZs4ft27dz+PBhRo0aVZRhi2xIMbEQQojSwqQoL96tWze6deuW7T5FUZg7dy6TJ0+md+/eAKxatQoXFxe2bt3KwIEDuXLlCjt37iQwMJDmzZsDsGDBArp3787//vc/KleuXJThi6dIMbEQQojSwmA1N2FhYURFRdGpUyfdc/b29rRs2ZLjx48DcPz4cRwcHHSJDUCnTp0wMjLi5MmTOV47NTWV+Ph4vU08H1l2QQghRGlhsOQmKioKABcXF73nXVxcdPuioqKoVKmS3n4TExMqVKigOyY7M2fOxN7eXrdVq1atkKMvXx4nPyY8LhyQBTOFEEKUfGXybqlJkyYRFxen2yIiIgwdUqmWOWrj6eCJvYW9gaMRQgghcmew5MbV1RWA6Ohoveejo6N1+1xdXbl//77e/oyMDGJiYnTHZMfc3Bw7Ozu9TRTc2aizADRxbWLgSIQQQohnM1hy4+XlhaurK/v27dM9Fx8fz8mTJ/H19QXA19eX2NhYgoKCdMfs378fjUZDy5Ytiz3m8upM5BkAmro1NXAkQgghxLMV6d1SCQkJhIaG6h6HhYVx7tw5KlSogLu7O+PHj+frr7+mZs2aeHl5MWXKFCpXrkyfPn0AqFOnDn5+fowcOZIffviB9PR0xowZw8CBA+VOqWIkIzdCCCFKkyJNbk6fPk379u11jydMmADA0KFDWblyJZ988gmJiYmMGjWK2NhY2rZty86dO7GwsNCds3btWsaMGUPHjh0xMjKiX79+zJ8/vyjDFk9JSk/i6sOrADRxk+RGCCFEyadSFEUxdBBFLT4+Hnt7e+Li4qT+Jp9O3jlJq2WtqGRdiaiPolCpVIYOSQghRDlR0M/vMnm3lCg8T09JSWIjhBCiNJDkRuTqbKTU2wghhChdJLkRudKN3Ei9jRBCiFJCkhuRowxNBhfvXwRk5EYIIUTpIcmNyNHVh1dJyUjB1swW7wrehg5HCCGEyBNJbkSOMuttGrk2wkglvypCCCFKB/nEEjk6dfcUAM3cmhk4EiGEECLvJLkROTp1T5vctKwiS10IIYQoPSS5EdlKzUjVrQbesqokN0IIIUoPSW5Ets5HnydNnUZFq4p4OXgZOhwhhBAizyS5EdnKrLd5ocoL0plYCCFEqSLJjciWLrmp/IKBIxFCCCHyR5Ibka2Td08C2pEbIYQQojSR5EZk8Tj5MdceXQMkuRFCCFH6SHIjssgctfF29MbJysnA0QghhBD5I8mNyOLw7cMAvOjxooEjEUIIIfJPkhuRxZHwIwC0c29n4EiEEEKI/JPkRuhJyUjR3SklIzdCCCFKI0luhJ5Td0+Rpk7DzcYNb0dZCVwIIUTpI8mN0PN0vY007xNCCFEaSXIj9Ei9jRBCiNJOkhuhk65O51jEMQDaeUhyI4QQonSS5EbonLhzgoS0BCpaVaRepXqGDkcIIYQoEEluhM6uG7sA6Fy9M0Yq+dUQQghROsknmNDJTG66enc1cCRCCCFEwUlyIwB4mPSQoHtBAHTx7mLgaIQQQoiCk+RGALDnxh4UFBq6NMTN1s3Q4QghhBAFJsmNAP6ZkupSXUZthBBClG6S3AjUGjU7QncA0LWG1NsIIYQo3SS5EQREBHA/8T6OFo685PGSocMRQgghnovBkxt/f39UKpXeVrt2bd3+lJQURo8ejZOTEzY2NvTr14/o6GgDRlz2bL6yGYBePr0wNTY1cDRCCCHE8zF4cgNQr149IiMjddvRo0d1+z788EP+/PNPNm7cyKFDh7h37x59+/Y1YLRli6IouuSmbx35uQohhCj9TAwdAICJiQmurq5Zno+Li2PZsmWsW7eODh06ALBixQrq1KnDiRMnaNWqVXGHWuYE3gskIj4CGzMbuQVcCCFEmVAiRm6uX79O5cqVqV69OoMHDyY8PByAoKAg0tPT6dSpk+7Y2rVr4+7uzvHjx3O8XmpqKvHx8XqbyN6my5sA6FGzBxYmFgaORgghhHh+Bk9uWrZsycqVK9m5cydLliwhLCyMF198kSdPnhAVFYWZmRkODg5657i4uBAVFZXjNWfOnIm9vb1uq1atWlG/jVJJrVGz9uJaAAbUHWDgaIQQQojCYfBpqW7duun+u2HDhrRs2RIPDw9+++03LC0tC3TNSZMmMWHCBN3j+Ph4SXCysefmHu4+uYuTpRM9fXoaOhwhhBCiUBh85ObfHBwcqFWrFqGhobi6upKWlkZsbKzeMdHR0dnW6GQyNzfHzs5ObxNZrTy3EoA3GryBmbGZYYMRQgghCkmJS24SEhK4ceMGbm5uNGvWDFNTU/bt26fbHxISQnh4OL6+vgaMsvR7nPyYrVe3AjCs8TDDBiOEEEIUIoNPS02cOJGePXvi4eHBvXv3+PLLLzE2NmbQoEHY29szYsQIJkyYQIUKFbCzs2Ps2LH4+vrKnVLPac2FNaSqU2no0pAmrk0MHY4QQghRaAye3Ny5c4dBgwbx6NEjnJ2dadu2LSdOnMDZ2RmA77//HiMjI/r160dqaipdu3Zl8eLFBo66dNIGCBkAACAASURBVNMoGuafmg/AqKajUKlUBo5ICCGEKDwqRVEUQwdR1OLj47G3tycuLk7qb4BtIdvovaE3DhYORHyo7XEjhBBClDQF/fwucTU3oujNPTEX0I7aSGIjhBCirJHkppwJuhfEgVsHMFYZM+aFMYYORwghhCh0ktyUM1MOTAG0t39Xs5feP0IIIcoeSW7KkYDwAHaE7sBYZcyXL31p6HCEEEKIIiHJTTmhKAqTD0wG4O0mb+NdwdvAEQkhhBBFQ5KbcmLzlc0cvHUQM2MzJrebbOhwhBBCiCIjyU05kJCWwPhd4wH4tM2nuNu7GzgiIYQQouhIclMOfHXwK+7E38HLwYtJbScZOhwhhBCiSElyU8YdvHWQOcfnALCg2wIsTQu20roQQghRWkhyU4Y9THrI4M2DUVAY0WQEPWr1MHRIQgghRJGT5KaMSlenM2jTIO49uYePkw/z/OYZOiQhhBCiWEhyUwYpisKo7aPYe3Mv1qbW/Nr/V6zNrA0dlhBCCFEsJLkpYxRFYcKuCaw8txJjlTG/DfiNRq6NDB2WEEIIUWxMDB2AKDxqjZr3/3qfn878BMCSHkvoXrO7gaMSQgghipckN2VETHIMgzYNYveN3ahQ8XOvn3m7yduGDksIIYQodpLclAHnos7R/7f+3Hh8A0sTS1a9uor+dfsbOiwhhBDCICS5KcXUGjX/O/Y/phyYQromHU8HT7a+vlVqbIQQQpRrktyUUlcfXmXUn6M4En4EgD61+7C051IqWlU0cGRCCCGEYUlyU8okpycz8+hMZh2dRbomHRszG+b7zWdY42GoVCpDhyeEEEIYnCQ3pYRG0bAxeCNf7P+CG49vANCjZg8Wdl+Ip4OnYYMTQgghShBJbko4RVHYFrKNKQemcPH+RQCq2FZhfrf5vFr7VRmtEUIIIf5FkpsSSlEUdt3YxZQDUzh97zQAduZ2TPSdyPhW47E1tzVwhEIIIUTJJMlNCXTw1kEm759MQEQAANam1oxvNZ6PfD/C0dLRwNEJIYQQJZskNyXIiTsnmLx/MvvC9gFgYWLB6Baj+bTNpzhbOxs4OiGEEKJ0kOSmBAi+H8wX+7/gj5A/ADA1MmVUs1F8/uLnVLatbODohBBCiNJFkhsDuh17G/9D/qw6vwqNosFIZcSwRsOY+tJUPBw8DB2eEEIIUSpJcmMAj5MfM+PIDOafmk+aOg2AfnX68XWH/2vvvsOjqNY/gH83nQAJLRAwMVSp0otB6UhQ1IsFUX8qKAgqFsqlKB0RUET0IgoWwIsoVUBFSqQLQaWEGjrSkyCQBAKpe35/vHcy2WQ3ZJPdzGbz/TzPPLM7OzP77uzs7DtnzjkzGfUq1TM4OiIiouKNyU0RSstMwxd/fYFJ2ybh2u1rAIDONTpjapepaH1Xa4OjIyIicg9MboqAUgo/xvyIkb+NzOqAr2FQQ3zU7SNE1IpgXzVEREQOxOTGiczKjB9jfsTkbZOxP24/AKBK6Sp4r9N7eKnZS/Dy4OYnIiJyNA+jA3BHl25cwge/f4D6s+uj17Je2B+3H2V8ymBs+7E48eYJvNLiFSY2znboENCpE+DjAzz0EHD5stERAbdvAydPAtHRQEaG/csnJQHXrwNKOT42Ikczm2VwlvR04O+/geRk571HSZWWBly4IONiqtj8w86ePRvTp09HbGwsmjRpglmzZqF1a2PrqZzZ/jNuXI/FrcwUxN26glMJZ/BX7B4cvhqDTBPgYQJa+JbBC0374MVmfVG+dEUg9hrgkQB4eMhgMukDkP/n1l7z8JA/c0/Pot0Q9sjMBOLigPPnZbh2TX5AOYfUVBlfvw6cOSOJQVAQ0LEj8PjjwD33WF9/SgowZQowbZoc/ABg3TrggQeADRuAWrVk2uHDwJIlMjabAS8vwNsbqFwZeOIJoF07fbsCEuelS0ClSkCVKvJaYqKs++efgaNHJb6uXYFevWS+XbuAbdtk2L9f1qGpWBEYMAAYNAi46y6ZdvMmcOqUvM+1a/LZr1wBDh4E9u2TAzkAVKsGdOgANG0KlC0LXL0q2+jMGdkHQkKAsDD5LF5e+lChAhAeLrEB8j7r10uciYmy71SoAFStKuvIyJB54uIAPz9ZrmpVoFw52bbp6bIdfH1lWqVKMvj4AKdPA0eOSJJ5+DCQkCDz1a4NtG4t30ft2rK8UsDFi8Dx47Lc6dPymdLSJJ7q1WWoXFm2x/HjQEwMcPasbPP27YFu3YAaNSz3hfR0Wc8//8iQkACULw/UrQsEB1vOm5Ym22HDBvkuU1KAu++WGKtWld/UpUvAsWPy3jduyD4YEQE8/LDElv19jx2T9w4MlG1386Zsh717ZZtkZkoc7dvLdxkSom+HI0dkP7hxQ5a7fVtiadkSaNVKvsvs+/vhw7J/7d8v265iRaBRIxn//TcQGwuUKSPx1q8v6woMBAIC9GOFUsCtW/IdeXnJNli4ENixA/D3Bx55BOjeXfavq1dlf6lUSdZVrpz+G9m5E9i4EfjtN/kcHh6yv951l7xHUpJs6zp15Lfcpo3+R3r+vKw3IEAGLUZAPx4kJMh69+8HDhyQ6Z6ewJNPAhMnAvVsNMi4elW2k8kk31XlyvIdbN8ObNoEnDghv7f0dNnvU1Nlu3t76/GEhckJU48esrxSEsuvvwJr1kg8SskxpmVLfQgLk+nXr8vveNs2ed+4OPmt1KsH3HuvPr+2v124IOvcs0f2mbQ0mT8sTPbLunVlCAnJfcxPTJT32rtXhj17ZN8oVUqOAa+9Jt+pR47yjcuXgcmTgW++kW3g5SUxtWsn+2rt2vK7O3JE9isfH9nXOnTI/ZsymEkp1z8NXLJkCV588UXMmTMHbdq0wSeffIJly5bh2LFjqJz9oGJDUlISAgMDkZiYiADtx+IAB2qVQePTLnjW4OkpB6nsg4+PjP395Q+xTBkZa4P2PPtBJedBpnRpOThkZsoPLTERiI+XPxxt/M8/ciDR/lSuXpUDWvYDtSM0awY0aSI/qKAg+dM6dw6YO1cvpXn0UeCtt4BXX5U/iypVgD595OC7Z0/e6+/YERg2TGJetAhYu1Y+NyDbIThYfuT2lsD4+8sB5eZNee7pKX/cSUmy/YpC3bryPR49WjTvZ0tQkPwxxsbK/lFYNWpIMnT9uux7SUm25w0LA+67T7b9gQPyZ6N9J/YymYDmzeW9L1+WxEZLrPMrOFjivXUr7/kCAuSPxNNTkrxjx/T9siBKl5bh2jV9X/bwsK/EJSBAji1Ftf9m5+Wlx+3pCTzzjCQKpUtLUnnhgiTwf/zhuBJPDw95j3/+kWTUaB4ekmxVqSK/p3Pn5ETnTho1AkaPlpLtK1eAL78EPvtMP0ZrJx751bixLN+uXcE+hw0F/f8uFslNmzZt0KpVK3z22WcAALPZjNDQULz55psYNWpUrvlTU1ORmpqa9TwpKQmhoaEOT27W3nMX6sfHw8Nsghc84A0v+Hl6wd/bG54KerHsnYbiQivJKOwuo53NhYbKH5yWfGlD9udly8qfVpky8oNdswaIjMz7gB4SAnzyiZTAmEzy59m9u5ztaby85Aysc2d5P+2M7eBB4LvvrBfHVqggZ47Zv7N69YDHHgPuv1/iW7lSzsyUkrPVDh1kuO8+/Yw5MxP45Rdg5kyZN7uKFWW+ChVk0EoamjWTkhofHzlQ79ghZ0+pqXJAq15dtpPJJAe3s2clsczMlCEjQz/j0phMUhLQpYscGFNT9RKqCxfkrLVqVfnjTUmRg/nly/In7O0tg1LyWkKCvH79uv456teXA2ijRrL+mzel1GPnTuDPPy23sacnULOmnPXWqiUHa19fSZpPn5bYr1yR/aVWLVl3zZryWSMjpdTFWqJpMkksFSvKdrp6VdZn7XdXqRLw4IPyXfn5yXuePCln2GazfIZ77pHvvHx5+Qxr1siZcU4BAbLdEhNlu/r6ynLa9+jtLWfWW7bI8lo8Xl5ydnzPPfIeZcrId37qFPD775alf9n3mSZNZLjnHr20LyFB9omqVSUOrdTp8mWJKS/e3lIa9uST8r0uX67/fipWlP04Pl62Z3a1aslvqmtXoG1b2f7nz0sS4OmplxZFR0vpztGj8vlCQ+V3W768JLpJSRJzUpKsQzsulC4tv4dGjYAWLeT9DhwAxo6VEtS81Kwpx54rV2TdANCggZTGtGghn8vHR74DHx8p5cjI0GM5cEBOcrJ/376+svzDD8vv3MtLfmO7d+uD9l4+PkDDhrJdOnaU2G/dkvmjo4G//pJtrP0utH2mZUvZZ8qUkcTjzBk9sT11ynYiHRYmSXfz5vL56teXWL7/Xk4CtbhyCg8H3n9fPs/Zs5L4a0NcnOzXjRrJd3bzpsQdHS3LHjokn9GB3Da5SUtLg7+/P5YvX46ePXtmTe/Tpw8SEhKwevXqXMtMmDABEydOzDXd0clNq1ay7+YUGCj/q/37y36S78ZQSulDXs9tvWY265d0tCHn81u39FKUGzesD9kPLNpgLZnw8JADQuXK8qcTFCR/ENqfifY4MNCydKhcOcuidXtduSJ/aH//LT+2K1fkTzUgQA4yzz4rB5LsEhMl4YmNlbOuXr0kXmvOn5ei2XXrJN5HHwVefFEODqmp8r4XL8qf0N13514+MVG2V/nyd/7ytUQkIEDWVd7J9w775x9JBDw9ZQfWLlE5ilakX7p03vOlpMgfcGqqJHG1a+f+zuyRlCQH2du3ZX3avleuXO4i+xs3ZN6oKNl/atSQP6jGjXMX0+fH+fOS6Ny6Jb+F+vXlwJ/fH76WeFSsKMvZ2g6ZmfInsn27/PGFhUlCU62aHQeZ/0lN1X/nycmyzcqVk+NFSor8ZsuUsVxGO9Zkf6/kZNmHb92SxKqsgTf03bFDks1LlySu27dlm7ZrJyc3ISH6vCkp8lkL8n9w9qz8iQcESNLg7297XqXkN2Ey6VUR8pKZKXGnpcn3caf5MzIkyYyLk+H6dfmc9evn/du+fl2Oh99+K5/H21suFQ8dKid99u5P8fHA1q3AU0/Zv+wduG1yc+nSJdx1113YuXMnwsPDs6aPGDECW7duxR9//JFrmaIqudm9W06CEhLk/27/fjkRy34y06QJMGaMJDsFOW66BO16fPazKG9v+VG7cv0eIiKyTTu2e3sX7uTCiQqa3BSbCsX28PX1ha+vr9Pfp2XL3NMyM+WE8JtvpL7q/v1SSNCwITBhgpTyFrtubUwm/do8ERG5B+3Y7oZcviyhUqVK8PT0RFxcnMX0uLg4BLtY7WxACjIeeACYP1+qLYwbJyW8hw9LktOqlTTIcO3yMiIiouLL5ZMbHx8ftGjRAhs3bsyaZjabsXHjRovLVK6oQgVpnfj335LklCkjjXQiIqQOp5UrakRERFRILp/cAMDQoUPx1Vdf4dtvv0VMTAxee+01JCcn46WXXjI6tHwpV06SnNOngcGD5dLm5s3SKOPxx6VUh4iIiByjWCQ3vXv3xkcffYRx48ahadOmiI6Oxrp161ClShWjQ7NLUJC0/j1+HHjpJalgvGqVNN7p00fvo42IiIgKzuVbSzmCszrxK6yYGGlJ9eOP8tzbW/qbGz1autQgIiIqyQr6/10sSm7cVf36wIoVUvemSxfpi2nWLOnbaexY230sERERkW1MblxA69bSWedvv0lrquRk6UOuZk3go48cd8cCIiKikoDJjQvRWlCtWCGlOteuAcOHyz3mvvzS/tvVEBERlURMblyMySS9GR84AMybJz3yX7wIDBwoHQGuWME+coiIiPLC5MZFeXlJi6rjx+UWIEFBwIkTcuuO7t3lMREREeXG5MbF+foCb78tfeSMHSvPN2yQm7KOG2f95tVEREQlGZObYqJMGWDSJLkZbffuktS89x7Qpo3cXJmIiIgEk5tipnZt4NdfgaVLgYoVgehouYHn9OmA2Wx0dERERMZjclMMmUxyE85Dh4BHH5VSnBEjpCJyUpLR0RERERmLyU0xFhwMrF4NzJ0r96tavVr6zDl1yujIiIiIjMPkppgzmYABA4Dt24GQEODYMaBdOynVISIiKomY3LiJ1q2Bv/6Sm3Bevgx06MC7jRMRUcnE5MaNBAcDW7ZIonPtmrSqunDB6KiIiIiKFpMbN1OhgrSmqltXEpsePYBbt4yOioiIqOgwuXFDFSsC69YBVarIbRwGDzY6IiIioqLD5MZNVa8OLFokFY6/+gpYvNjoiIiIiIoGkxs31qULMGaMPH79dSA+3th4iIiIigKTGzc3bhzQtClw/Trw738bHQ0REZHzMblxc15e0smfyQQsXAhs2mR0RERERM7F5KYEaN0aeO01eTxsGO9BRURE7o3JTQkxcSJQtqzcaHPZMqOjISIich4mNyVEpUp6nZuxY4H0dGPjISIichYmNyXIkCGS5Jw4Ic3EiYiI3BGTmxKkbFlg+HB5PGMGoJSx8RARETkDk5sSZsAAoEwZuWv4hg1GR0NEROR4TG5KmHLlgH795PGMGcbGQkRE5AxMbkqgt98GPDyAyEjg8GGjoyEiInIsJjclUI0awGOPyeNvvjE2FiIiIkdjclNCaZemFi4E0tKMjYWIiMiRmNyUUN27A1WrAv/8A/z8s9HREBEROY6hyU316tVhMpkshmnTplnMc+DAAbRr1w5+fn4IDQ3Fhx9+aFC07sXLC+jTRx7z0hQREbkTw0tuJk2ahMuXL2cNb775ZtZrSUlJ6NatG8LCwrBnzx5Mnz4dEyZMwJdffmlgxO7j5ZdlvH49EB9vbCxERESO4mV0AGXLlkVwcLDV1xYtWoS0tDTMmzcPPj4+aNiwIaKjo/Hxxx9jwIABNteZmpqK1NTUrOdJSUkOj9sd1KkDtGwJ7N4NrFolfeAQEREVd4aX3EybNg0VK1ZEs2bNMH36dGRkZGS9FhUVhfbt28PHxydrWkREBI4dO4br16/bXOfUqVMRGBiYNYSGhjr1MxRnTz0lY95Mk4iI3IWhyc1bb72FxYsXY/PmzRg4cCCmTJmCESNGZL0eGxuLKlWqWCyjPY+NjbW53nfeeQeJiYlZw/nz553zAdzAk0/KePNmqVxMRERU3Dk8uRk1alSuSsI5h6NHjwIAhg4dio4dO6Jx48Z49dVXMWPGDMyaNcviklJB+Pr6IiAgwGIg62rXBpo2BTIzgdWrjY6GiIio8Bxe52bYsGHo27dvnvPUrFnT6vQ2bdogIyMDf//9N+rWrYvg4GDExcVZzKM9t1VPh+zXqxcQHQ0sX673f0NERFRcOTy5CQoKQlBQUIGWjY6OhoeHBypXrgwACA8Px+jRo5Geng5vb28AQGRkJOrWrYvy5cs7LOaSrmdPYPRouTR16xbg7290RERERAVnWJ2bqKgofPLJJ9i/fz9Onz6NRYsWYciQIXj++eezEpfnnnsOPj4+6NevHw4fPowlS5bg008/xdChQ40K2y3Vrw/cfTeQmgps3Wp0NERERIVjWHLj6+uLxYsXo0OHDmjYsCHef/99DBkyxKIPm8DAQGzYsAFnzpxBixYtMGzYMIwbNy7PZuBkP5NJeiwGgHXrjI2FiIiosExKKWV0EM6WlJSEwMBAJCYmsnKxDatWAY8/Ln3fHD9udDREREQF//82vJ8bcg2dO8stGU6cAE6dMjoaIiKigmNyQwCAgADg/vvl8fr1xsZCRERUGExuKEtEhIwjI42Ng4iIqDCY3FCWTp1kvH07YDYbGwsREVFBMbmhLC1aSB83V68CMTFGR0NERFQwTG4oi7c30LatPGZ/N0REVFwxuSEL7dvLeNs2Y+MgIiIqKCY3ZKFDBxlv3Qq4fw9IRETkjpjckIXWrQFfXyA2Fjh50uhoiIiI7Mfkhiz4+QFt2shjXpoiIqLiiMkN5aJ15rdrl7FxEBERFQSTG8pFK7n54w9j4yAiIioIJjeUi5bcHD4M3LhhbCxERET2YnJDuQQHA3ffLb0U795tdDRERET2YXJDVvHSFBERFVdMbsgqJjdERFRcMbkhq7InN+zMj4iIihMmN2RV8+aApydw+TJw4YLR0RAREeUfkxuyyt8fuPdeecxKxUREVJwwuSGbmjeX8b59xsZBRERkDyY3ZJOW3Ozda2wcRERE9mByQzY1ayZjltwQEVFxwuSGbGrSBDCZgEuXgLg4o6MhIiLKHyY3ZFPp0kDduvKYpTdERFRcMLmhPLHeDRERFTdMbihPWr0bJjdERFRcMLmhPLE5OBERFTdMbihPTZvK+PRpICHB2FiIiIjyg8kN5alCBaB6dXkcHW1oKERERPnC5IbuiPVuiIioOHFacvP++++jbdu28Pf3R7ly5azOc+7cOfTo0QP+/v6oXLkyhg8fjoyMDIt5tmzZgubNm8PX1xe1a9fGggULnBUy2aAlNyy5ISKi4sBpyU1aWhp69eqF1157zerrmZmZ6NGjB9LS0rBz5058++23WLBgAcaNG5c1z5kzZ9CjRw906tQJ0dHRGDx4MPr374/169c7K2yyQqt3s3+/sXEQERHlh0kppZz5BgsWLMDgwYORkKM26tq1a/HII4/g0qVLqFKlCgBgzpw5GDlyJK5cuQIfHx+MHDkSa9aswaFDh7KWe+aZZ5CQkIB169blO4akpCQEBgYiMTERAQEBjvlgJcj588DddwNeXsDNm4Cvr9ERERFRSVDQ/2/D6txERUXh3nvvzUpsACAiIgJJSUk4fPhw1jxdu3a1WC4iIgJRUVF5rjs1NRVJSUkWAxVcSIhULM7IAI4cMToaIiKivBmW3MTGxlokNgCynsfGxuY5T1JSEm7fvm1z3VOnTkVgYGDWEBoa6uDoSxaTSe4zBbDeDRERuT67kptRo0bBZDLlORw9etRZsebbO++8g8TExKzh/PnzRodU7Gn1bpjcEBGRq/OyZ+Zhw4ahb9++ec5Ts2bNfK0rODgYf/75p8W0uP/dejo4ODhrHJfjdtRxcXEICAhAqVKlbK7b19cXvqwY4lBMboiIqLiwK7kJCgpCUFCQQ944PDwc77//PuLj41G5cmUAQGRkJAICAtCgQYOseX799VeL5SIjIxEeHu6QGCj/sreYUkouVREREbkip9W5OXfuHKKjo3Hu3DlkZmYiOjoa0dHRuHnzJgCgW7duaNCgAV544QXs378f69evx5gxYzBo0KCsUpdXX30Vp0+fxogRI3D06FF8/vnnWLp0KYYMGeKssMmGevUAHx8gMRE4e9boaIiIiGxzWnIzbtw4NGvWDOPHj8fNmzfRrFkzNGvWDLt37wYAeHp64pdffoGnpyfCw8Px/PPP48UXX8SkSZOy1lGjRg2sWbMGkZGRaNKkCWbMmIGvv/4aERERzgqbbPDxARo2lMe8NEVERK7M6f3cuAL2c+MYL78MzJ8PjB8PTJhgdDREROTuil0/N1T8sDk4EREVB0xuKN/YYoqIiIoDJjeUb1rJzdmzwPXrxsZCRERkC5Mbyrdy5YDq1eXxgQOGhkJERGQTkxuyCy9NERGRq2NyQ3ZhckNERK6OyQ3ZhckNERG5OiY3ZBetUvHhw0BamrGxEBERWcPkhuwSFgYEBgLp6UBMjNHREBER5cbkhuxiMvHSFBERuTYmN2S37HcIJyIicjVMbshuLLkhIiJXxuSG7JY9uXH/264SEVFxw+SG7Fa/PuDlJbdgOH/e6GiIiIgsMbkhu/n6Ag0ayGNemiIiIlfD5IYKhPVuiIjIVTG5oQJhiykiouInJQWIjzc6CudjckMFwpIbIqKit38/8OKLwDPPADt22Lfs8eNArVpAcDAwebJz4nMVJqXcv71LUlISAgMDkZiYiICAAKPDcQvXrgEVK8rjhATptZiIiJznzz+BDh2k9EUzezbw+uv5W757d2D9enns6SnJTs2ajo/TkQr6/82SGyqQChWAu++Wx7w0RUZKTWWXBOT+btwAeveWxKZzZ+C552T6G28AkZF3Xv7MGUlsPDyA2rWBzExg4ULnxmwkJjdUYM2by/ivv4yNg0qmjAygTx/Azw9o0QK4dMnoiIicZ/p04O+/gerVgR9/BL77DujXTxL7gQOBW7fyXn7pUhl37AiMGSOPV692YsAGY3JDBRYeLuOoKGPjoJLpq6+A//5XHu/bB7zyirHxEBVUVBSwYAEQF2f99dhYYMYMeTxjhlQDMJmAmTOBkBAplZk7N+/3WLJExs88Azz4oDzev19KhNwRkxsqMC252bmTlwWoaN24AUyYII8HDJBOJX/9Fdi2zdCwiOw2aRLQti3w0ktAw4bA9u2555k6VUpm7rsPePxxfXrZssC4cfJ4xgy5RGvN8eNyAuDlBTzxBFCtGhAWBpjNwB9/OP4zuQImN1RgLVvKj+XyZeDcOaOjoZJk+nRpzlqnDvDZZ0D//jL9vfeMjYvIHrt3A+PHy+PgYODqVeCxx4BTp/R5LlzQS2Xee09KbLJ78UVJVi5eBJYvt/4+WqnNgw/qDUHatpXxzp2O+SyuhskNFVipUkCzZvLYXX8g5HgZGYVb/tIlvYh+2jTA2xsYOVJaf/z2m5yhEhUHw4fL+PnngdOnpWQmIUFKV5KT5bX335cSmfbtgS5dcq/D1xd47TV5PHu29fdZvFjGvXvr07Tkxl2rFTC5oUJx9+yfHCchQZqienvLGeT16/lbLiYG2LoVuHlTLn++8YYU0d9/v15EX726fuCePt0p4RM5VEwMsGWLJOVTpsjJ4vLlQJUqwIEDUlk+MlIvtZk0KXepjaZ/f/ldRUXlTu4PHACOHAF8fICePfXpLVvKeO9e96xWwOSGCoWViim/+vfX+9j47Tdp4ZHXQVVLZBo0kBYeFSoA9eoBK1fKgfyzzywP9tpZ8NKlchZMJUdhSwON8O23Mn74YSA0VB7fdZckOF5ewIoVQLdu8jvo21f6t7ElOBh48kl5/Pnnlq99/bWMH33Usj+yxo2lWXh8vFRYLoiYGGlOfvlywZZ3JiY3RnGpFAAAIABJREFUVChayU10tF6MSpTTn3/KwdrTE/j0Uzl4L1sG/Pyz7WW+/14vZq9WDUhPl4qRgFyW0nrJ1jRtCkRESP8dr7wiY3J/8+ZJxdouXYrPMSh7HzN9+li+9sAD8ruoUEGed+oE/Oc/d17noEEyXrRILxW9eVN/nwEDLOf39wfq1pXHd7qUm56e+/f03/8C994rdX4aNZJEx6WoEiAxMVEBUImJiUaH4pZCQ5UClNqwwehIyFX17Cn7SJ8+8nzkSHler55S6em55z97VqnAQJln4kSlzGaljh9Xat48paKibL/P8eNK+fvLcg8/rNSZM874NGLhQqWeflqpRYuc9x6Ut9hYpUqVku8bUGrSJKMjyp916yTeChWUSkmxPk9ysuy/ZnP+1mk2K3XvvbLeceNk2tix8rx2baUyM3Mv8+yz8vrkybbX+8EHso3LllVq9GilbtxQ6pdflPL0lGW131ubNvmP1R4F/f9mckOF1qeP7NwjRhgdCbmiQ4dk/zCZlDpyRKYlJChVqZJMnzPHcv7MTKU6dpTX7rvPevKTl1WrlPL11f/wWrRQqksXperXVyo4WKmBA5W6datwn+mHH/T1A0p9/33h1kcFM3u25fdQq5Zz/mDtcfy4Uv/3f0o9+aRSO3ZYn+eZZyTeN95w7HsvXSrr9faWRE/7HaxYYX3+Dz+U15980vrrK1ZYbl9ATjq0xOaFF5S6eFGp0qXzfp/CcLnkZvLkySo8PFyVKlVKBQYGWn9zINfwww8/WMyzefNm1axZM+Xj46Nq1aql5s+fb3csTG6c67//1f9EqGQ5d04OqOvX2z4Dff556wfQ//xHplepImeDmunTZXrp0kqdOFGwuPbtU6pr19wHZm148cWCrVcppTIylKpTR9YTHCzjkBCl0tIKvk4qmIceku0/dqxegrN7t3Hx/POPUtWq6fuZp6dS33xjOc+1a3rS4ehYzWa9lFQbHn/cdsIXGSnz1KxpfV3168vrQ4ZI4lKzpr7eJ57Q9/lJk5QaNEhK0hzN5ZKbcePGqY8//lgNHTo0z+Rm/vz56vLly1nD7du3s14/ffq08vf3V0OHDlVHjhxRs2bNUp6enmrdunV2xcLkxrkuXtTPzK9eNToacgazWakff5QD2BNPSMlKjRq5z5r/+styuZMn9bO8nAfy1FRZBlBq/HiZtm2bnHUCSn35ZeHjvnBBkq+FC5X67Te5hGQyyfr37CnYOn/8UZYvV06pK1ckOQMkyaeic/OmniQcOqRUr17y+J13jItp2DCJoU4dSeYBpTw8lFq2TJ9nxgyZfu+9zillSk6WONq2laTP1kmHUpKMab/fhATL137/Xb/spP113r4tv6GffrJ+mcsZXC650cyfPz/P5GblypU2lx0xYoRq2LChxbTevXuriIgIu2JgcuN89eo5r1iSjHX7ttRfsVYC4uGhVLNm+iWmcuWk1ETzf/8n0x96yPq6tWJ0Hx8poi9TRp736uW8ywvaJYFXXy3Y8uHhln+i778vz1u2dFyMdGerV8t2r15d9pVvv5XnrVsbE09qqv47+PlniWnAAH3/3rhRSjq0k4KvvjImzpy0OpNbt1pO79dPpr/0kjFxaYptclOtWjVVsWJF1apVK/XNN98oc7YjWrt27dTbb79tscy8efNUQEBAnu+ZkpKiEhMTs4bz588zuXGyN9+UH0K/fkZHQo5kNivVt698t35+koDMnq3U4sVysNZ+UomJSt1/v8xXtapUCN6wQS8l2bvX9vqfesoyYerUSc4+nUUrii9fXv6Q7KGdzfr4KHX5skyLj9dLm7Indndy8qTUySiqM2B388orlvVWzp3TE+6cpRCOsGuXUtOmKRUTY/31Vav0S5VaPbGMDH3/LltWqX/9Sx5XquTcfdwejz4qMX3yiT4tM1MvkfztN+NiU6qYJjeTJk1Sv//+u9q7d6+aNm2a8vX1VZ9++mnW63Xq1FFTpkyxWGbNmjUKgLqVR43A8ePHW63Pw+TGedav1+tP8GDtPj79VP/DuFNruOvXlWrYUE9wAgLk8cCBeS+XkSGXoPr2Verzz+2vQGyvjAy9XsTPP9u3rPZH0L+/5XTtksjrr+dvPYsWyTYFZJ0ZGfbFUdKZzbKPAdLySKPVhfrpJ8e+39q1+uXVsmX1ivHZPf64vP7vf1tOT0lRqnNnywT+u+8cG19hjBsnMWktGZWSkxGt3ltel7WKQpEkNyNHjrSaNGQfYnKktXklNzmNHTtWhYSEZD0vaHLDkpuil5qq/5nt3Gl0NJTT0aNK/fqrXGK6k4wMOXh/8IFe8jJjRv7e59w5aXaqHcTbti18yyRneOMNie+VV/K/zPbteqJ39Kjla1ppUGDgnc/Is7cu0QZH1C8qSfbs0f98s+/TAwfK9MGDHfdet24pFRZm+X3lvMz6zz966d3Bg7nXcfu2XL588kmllixxXGyOsHKlxN2kiT5Nu9T62GPGxaUpaHJjVyd+w4YNQ0xMTJ5DzZo17e9s53/atGmDCxcuIPV/tzYNDg5GXI57wMfFxSEgIAClSpWyuR5fX18EBARYDORcPj7AQw/J459+MjYWsvT550D9+tITaps2ckdtW5YulVsZNGgg92tSCnj9dWDIkPy9V2iodOi4aJF00rd5s3Qr72oefVTGv/wid0bOS2IisGuXdAwISE/LWudnms6dgRo1ZN7vvst7fR98IJ3N3Xeffo+sd9+1/b1cugRMnAjMmmX7rs8lzerVMn7wQcDPT5+u3Xtp40bHvdcXXwBnzwIhIcD+/dKr79q1clsDzQ8/SEd3zZpJh3Y5+fnJd7x8OfD0046LzRG0+wMePqzvX+vWybh7d2NicggnJVtZ7Cm5mTx5sipfvnzW8xEjRqhGjRpZzPPss8+yQrGL+v57yfbvucf4viZIvoMPPshdCThnsbkme98t/v7Sx8zMme55mTElRa+8/Oef1ufZtUsqCWffdnfdpVRcnPX5P/5Y5qlcWZr7WnPpktRd0uoypKXpl1KmTs09/9Gjsj7t/Z9+umCf152YzXoru4ULLV+Lj9e3la3vyR6Zmfp7zZ0r055+Wp4//7w+X4sWueutFBdms9Q/w/9aEF67pl+CO33a6OhcsM7N2bNn1b59+9TEiRNVmTJl1L59+9S+ffvUjf91aPHTTz+pr776Sh08eFCdOHFCff7558rf31+N07pWVHpT8OHDh6uYmBg1e/ZsNgV3YYmJel8Tf/xhdDQlW2ys3rkiID2L/vyzPC5TRurHZHfunH6pZODA/F2+Ku60prpjx+Z+LSpKv8ygJSw9e+Z9sE9N1VsNhoZKB2dz5lhepho0SL9cp50ALFwo0ypUUCopSZ/39GlJprT1aX84W7Y45vMXVzt26Jekbt7M/XrjxvK6Iy7/aHUJAwL09/rrL5nm6SkV56Oi9Erm8fGFf08jdOqkXx797jt5nKOhsmFcLrnp06eP1To5mzdvVkoptXbtWtW0aVNVpkwZVbp0adWkSRM1Z84clZnjNHHz5s2qadOmysfHR9WsWZOd+Lm4556zbMFARW/BAj1RMZmU+ugjmW426xV+P/zQchmtFcf997tnSY01CxbIZ27a1HJ6crLeXPehh6SOTH5FR+sVXbUhJET+JPft0xOUjRv1ZdLTpbQTUOq112Ta3r2yHCAdqcXH682Kn3uu8J+9OIuIkO3Qt6/1199+23JbFobWId6bb1pO15KBIUP0JNlWPMWBdpuGXr301l2jRxsdlXC55MaVMLkpOto9UypWtL+Z7Z0sXSqXSkpKyUJBaKUAgPT3sW2b5etffy2v1aiht9DRKhR6eUlnaCVFfLzeYuncOX261iV9aGjBmhTfuCGtdSZMsKyIqnU4Z62r+w0b9PmaN9eToPr19eTqjz9kWqlSzmnq7CquXZPSkAMHlDp/XkpMzGYZ5s2TbeDtLU3prdH25/r1CxfHuXP6/nH4sOVra9daJrAmk1L79xfu/Yyk7VvZh5wdchqFyU0emNwUnfR0/cw1x500CmXnTv1A4+jWEO7iwgW9HsmQIdZLYG7dkssfgPwBJyToJQRG9uxqFK1vns8+k+eJifr2KUAhcS43b+p9QAHS4eGVK9bnnTLF8s/lySelFY7GbNYve9kb2+bNSg0fLvXiXLVkLjVVqaFD9cQu++DlJU2wtefvvmt7PVev6q38CnM7AK00o2PH3K+ZzZYdW1q7tFmcmM1KNWpkmWC7Sr1JJjd5YHJTtCZOlB9ImzaOW6f2J6TVQfD2ds59TIqz11+XbXPffXn3mzJ8uF7vQ7uMWKuW63QqVpS0SsANGsif/rvvyvO6dR3b386JE9IB4J3WuX+/3MbBWj8qSklpEKBUjx75f+9ZsywThW7dLO/l5QrS0/XLO9rvvHJlSWqyx+7rK/vvnfoFatKkcPVu0tL0+4bZWkdamtyKY9Mm10kECmP7dqWCguRz26pkbwQmN3lgclO04uKkch3gmD5vsvcKe+mSUq1aFd+WCc5y7py+zTdtynves2f1it9akfr27UUTp6tJSND7Z/r3v/VtmMddYQx1+LCe3OesFG7Npk16iedDD0krOEA6lXN2Z4m2mM1KTZ4sJQUvvSQVhLUWSD4+kjBkn/fmTbk8dfBg/j6zUlKyCxT8FhvabUGqVHH85XVXlpbmejeAZXKTByY3Re+ll+w/w7RFq+yq9QqrnYk2a5a/5XfvlhszHjtW+FhcldZ5WYcO+Zt/yRK5hOXvn/uuxSWNdvlBGx56yLXPxBs0kDi//db662azJAObNulNfF98Uab/8Yd+6dKoSyna3eBzDt7ecgsDR9BuhVCvXsGW1yoMjxnjmHio4Jjc5IHJTdE7dky/dr5jR8HXExOjly5ovcJeuZL/fhhSUvSu9hs0cM9u7o8f14vvc1YgzktKSsk6K7UlLU1auvj5ySWb7PVcXNH48fJdP/JI7tcuX87dN0+bNpa9RGv9GXl7F30/Jikp+uWe116T1jne3lLBPXsLssLKXu9GuwdYfh05Ist5eFhWNCdjMLnJA5MbY/TvLweJdu0Kfib8/POyjn/9y3J6+/YyfdasvJfX+nbRBnv+/IuDmzf1PzNbd94m93LokH4JJ3vlZLNZqS5d9H29cmX5/VhL1rp2lXleeKHo4lZKKkJrdWq0yx+Zmc4pKWvaVN5r8WL7lnvrLVnOFW49QEV0+wUie4wbJ92Ob98u3frb69Ah6cYfAMaOtXwte/f5efn9d8vn69fbH4cr2r4deOop6RJ+926gYkXpJp7cX8OGQMuWQFoa8J//6NOXLZPbDvj5AUePAnFxwMKFsm/kNHWqjBctAk6dKpq4lQI+/lgev/UW4O0tjz08AJPJ8e/XsaOMt2zJ/zJJScC338rj1193dERUlJjckNOEhgLvvCOPhw0Dbt7M/7KJicBzz8kBsWdPoEULy9cfeUTGmzfnvd6//pJxu3Yy1u6ZUpwtXQp06ACsWAEkJABhYZLkhYUZHRkVlZEjZfzpp3Lfo6tX9ft/jRqV+95XObVsKfeCM5uBadOcG6smMhI4eBAoUwYYMMD571eQ5OaLL+TYU7++3LeKijEnlSS5FF6WMs7t20rVrCnFvCNH3nn+s2elpUNgoF60fulS7vmy318mr5YtWh8uK1boxfW2+hnZvl06bmvcWOr6uKKkJL2SaO/ecv8jo1q9kHEyMqQuDSD76wMPyOM6dfJ/F3btNgbe3kVTt+TBB+X93n7b+e+llHQGqLUUy8/v+dYtaR0FSO/V5BpY5yYPTG6M9dNP+kF03z7r81y6JLds0Jriagfq3bttr1frZr1fP+uvJyfr6/rnH/3WA8uX55733Dnp40Gbv3Nn+z9nUfjkE33bMKkp2U6fVqpSJcs+YPbutW8dWqugnLcXcLStW+V9PD2LthLzI4/I+w4ffud5Z8+Wee++2/WaQ5dkrHNDLuuRR+TSUnq6XGq6dUt/7coV4N//BmrWBD77TOoRdOoErF0r9QZyXo7KrkcPGa9fL4f3nE6elHH58lLvQCum3rzZcr6UFOCJJySWgACZtmkT8PffBfm0zqMUMHu2PB42DPDyMjYeMlaNGlLf6q23gJdeAnbsAJo1s28dY8bIeO5c4MABx8cIyH6rXUYbMEDiLir9+8v466/zrluUkQFMny6Phw/X6wNRMeakZMulsOTGePHxehPQ7t3lcsqIEfoNHrUec+1pDpqcrJf0nDiR+/Vlyyx7Sl6+XG8Snt2rr8r0ChWUOnNG+orJT0usoqY1UfXxcb0eZql4MpulLypA7sXkjP3qxx9l/f7+9jfLLqz0dL23YkBaUFm7LYx2T7agoJLZU7crY8kNubSgIGDxYqBUKanUe999wIcfAsnJUjqzdq20bOrcOf/r9PeX9QC5S2MA4MQJGd9zj4w7dpSWGUeO6KU669cDc+bI4x9+AKpXB7p0kedRUfZ+Suf6+WcZd+4slTKJCstkAubPB6pWBWJipKQjeymo2Qy8956Ups6fb//6zWZpNQlIhefgYMfEnV9eXsDKlcD998vz6Gjg2WelMn72GLXWY0OGyHGFij8mN1RkOnSQlgv33y/JTqdOwE8/SYum7t0L1hy0UycZW0tujh+XsZbcVKyoJy5LlgDXrwP9+snzt94CunWTx23ayHjXLvvjcaaffpLxY48ZGwe5l6AgaUbu5SW/C637BQD4/HNJTrZsAV5+WU5I7LF6tXTpEBAgl5+NUKOGnDjFxgKvvCLT+vYFzp2Txz//LCc8AQFs/u1WnFSS5FJ4Wcp9bdkixcnBwbk7AmvbNncnXt98I9Nq15ZOurTKudmLoq9f14ux4+MdH7PZLHeAbt1aqXnz8rdMfLze4yp7TSVnmDxZv5/S9ety3y2twrJWGR9QasAA6Rzz5s07r1PrLDCvu3gXpfR0/Sa8Dz8sFYfr15fno0YZHR1Zw9ZSeWBy475SUqTbfGvNPbXWT9lbkCQmKlWxouX9bP74I/d669WT13/5xfExazfl04b83J5iwQL77qdFZK+UFLkbOiAtF8eM0e/PlJ6u1MSJlvtttWpyF3NbYmP1ptgnTxbd57iTmBi9rp5WH6dSpfzflJOKFuvcUInk66vXu9m+XZ9+/bq0fgKA2rX16QEBwLx5QOnSsuyXXwKtW+derzZtzx7HxzxrluXzoUOtt/bKbvVqGWs9MxM5mq+v3hrvs8+AyZPl8ZQpcslq3Di5hPP889JB56VL0goyJcX6+pYtk/osrVsDtWoVzWfIj3r1gA8+kMf798t45kygXDnjYiLHY3JDxZ7W+/C2bfq0I0dkHBoKlC1rOf9jjwHx8dKra9++1tfZvLmM9+51aKhISgJ27pTHUVFSefGPP4ANG2wvk5ys96z8+OOOjYcouy5dgGee0Z8/8IAkMJpHHpFbOhw4ILf+OHNGEiFrfvhBxtnX5yrefhv4/nv5/a9YIQkbuRcmN1TstW8v4+wlN4cPy7hhQ+vL+PtL6Y0tzkpuNm8GMjOlkvN99+nd0OfVBf769cDt21IxskkTx8ZDlNOXX0r/N8OGAatWWa/oX66ctKICpKVRcrLl62fPShJvMgG9ezs/ZnuZTNJqav586eOK3A+TGyr2wsOl2PzsWb0FhFZy06BBwdbZtKmMz5/XL285glZCo923Ztgw6TBsyxbbrbO0G/k98YRzbjBIlF3ZspK4fPSR9Ztual54QS43Xbsml3qzW7JExh06ANWqOS9WIluY3FCxV7q0XtKild7cqeTmTsqW1ZuQO7L0RktutGbnISF6kfiECbnr3pw6pfdvozVjJXIFnp568+4ZM6QHckD2Ya05uStekqKSgckNuQXt0pRW70ZLbgpacgM4/tLU6dPSeaCXl34rCEC6pvf2lstPY8YAqakyXeu2XinpB+hOd3omKmp9+gCVK0up6dKlMu2vv6ROjq8v0KuXsfFRycXkhtyCVql461Y50F6+LGeWjRoVfJ1acvPXX4WPDwAiI2UcHq7fwwqQpGXGDHk8ZYrUZ2jeXFp1rFghn+P99x0TA5EjlSolHWACwOjRwI0beourXr2AChWMi41KNiY35BbatZMSkWPH9NsptGhRuNsUPPCAjLdskUrAhZXzklR2b7wh9RaqVZOmtfv2SQ/Lvr7AggV6okXkat5+W25bcvas/GYWLpTpb7xhaFhUwjG5IbdQvrxeSVdrefTww4VbZ6tWQGCg9Jmze3fh1pWRIXcaB/Q4szOZ5M7OFy5Igvbrr9JS5eJFNlMl11amDPDf/0oifuCAXEZ95RX9NiZERvAyOgAiR/n3v+UGnIAcaG31YZNfXl5A165yaWjt2sIdrHfuBBISpPVJy5a25zOZpCKzVpmZqDho1w7480/gu++kBZV2zzYio7DkhtxG587A3LmSkCxbBoSFFX6dWo/A8+fb7ok1P375RcYPPyx1aIjcTePGcmPNgQPlxIDISExuyK0MGCAVdx11m4KnnwaCg6X/nKeekh6G7aWUfkfvRx5xTFxERGQbkxuiPJQqJRV6/fyANWukpVN8vH3r2LdP6tH4+QEREU4Jk4iIsmFyQ3QHERHSf061atLz8euv27e81nrkscekgjIRETmX05Kbv//+G/369UONGjVQqlQp1KpVC+PHj0daWprFfAcOHEC7du3g5+eH0NBQfPjhh7nWtWzZMtSrVw9+fn6499578euvvzorbCKrWrWSSsUmk1Qwjo62fD0lBXjySUleBg+WuyEDwK1benLzwgtFGzMRUUnltOTm6NGjMJvNmDt3Lg4fPoyZM2dizpw5ePfdd7PmSUpKQrdu3RAWFoY9e/Zg+vTpmDBhAr788suseXbu3Ilnn30W/fr1w759+9CzZ0/07NkThw4dclboRFY1bqx3Jz9liuVr48cDP/4odXI+/VQGQLqhv3pV+gF56KEiDZeIqMQyKZXzbjbOM336dHzxxRc4ffo0AOCLL77A6NGjERsbCx8fHwDAqFGjsGrVKhw9ehQA0Lt3byQnJ+MXrbkJgPvuuw9NmzbFHK23tjtISkpCYGAgEhMTEZC9a1giOx04IHfm9vSUTsvuuguIjQVq1pQ7d7dvL5ew/Pyk+XfPnlIZ+aOP5CaZRESUfwX9/y7SOjeJiYmokK0/7qioKLRv3z4rsQGAiIgIHDt2DNevX8+ap2vXrhbriYiIQFRUlM33SU1NRVJSksVA5AiNG0sCk5mp94T84YeS2LRpA2zeLD0Qp6RIr8LnzgGhocCrrxobNxFRSVJkyc3Jkycxa9YsDBw4MGtabGwsqlSpYjGf9jw2NjbPebTXrZk6dSoCAwOzhtDQUEd9DCK8+aaM586VkpzPP5fnEyYAHh7AN9/I/aEAqaPz2Wdy53IiIioadic3o0aNgslkynPQLilpLl68iO7du6NXr1545ZVXHBa8Le+88w4SExOzhvPnzzv9Pank6NlTLkNduSKXqFJT5ZYKWjPvkBBg1y652eW2bdJKioiIio7d/UgOGzYMfe/Qr33NmjWzHl+6dAmdOnVC27ZtLSoKA0BwcDDi4uIspmnPg4OD85xHe90aX19f+Pr63vGzEBWElxfwxRdSQdhsBsqWBWbNklIaTd26QLa680REVITsTm6CgoIQFBSUr3kvXryITp06oUWLFpg/fz48PCwLisLDwzF69Gikp6fD29sbABAZGYm6deuifPnyWfNs3LgRgwcPzlouMjIS4eHh9oZO5DDdukmF4W3bpCSnTh2jIyIiIo3TWktdvHgRHTt2RFhYGL799lt4ZruhjlbqkpiYiLp166Jbt24YOXIkDh06hJdffhkzZ87EgAEDAEhT8A4dOmDatGno0aMHFi9ejClTpmDv3r1o1KhRvmJhaykiIqLip6D/3067vVlkZCROnjyJkydPIiQkxOI1LZ8KDAzEhg0bMGjQILRo0QKVKlXCuHHjshIbAGjbti2+//57jBkzBu+++y7q1KmDVatW5TuxISIiopKlSPu5MQpLboiIiIqfYtHPDREREZGzMbkhIiIit8LkhoiIiNwKkxsiIiJyK0xuiIiIyK0wuSEiIiK3wuSGiIiI3AqTGyIiInIrTG6IiIjIrTC5ISIiIrfC5IaIiIjcitNunOlKtNtnJSUlGRwJERER5Zf2v23vbTBLRHJz48YNAEBoaKjBkRAREZG9bty4gcDAwHzPXyLuCm42m3Hp0iWULVsWJpPJYetNSkpCaGgozp8/z7uN3wG3Vf5xW+Uft1X+cVvZh9sr/5y5rZRSuHHjBqpVqwYPj/zXpCkRJTceHh4ICQlx2voDAgK48+cTt1X+cVvlH7dV/nFb2YfbK/+cta3sKbHRsEIxERERuRUmN0RERORWPCdMmDDB6CCKM09PT3Ts2BFeXiXiCl+hcFvlH7dV/nFb5R+3lX24vfLP1bZViahQTERERCUHL0sRERGRW2FyQ0RERG6FyQ0RERG5FSY3RERE5FaY3BAREZFbYXJzB7Nnz0b16tXh5+eHNm3a4M8//8xz/mXLlqFevXrw8/PDvffei19//bWIIjWePdtqwYIFMJlMFoOfn18RRmucbdu24dFHH0W1atVgMpmwatWqOy6zZcsWNG/eHL6+vqhduzYWLFjg/EBdgL3basuWLbn2K5PJhNjY2CKK2BhTp05Fq1atULZsWVSuXBk9e/bEsWPH7rhcST1eFWR7ldRj1hdffIHGjRtn9T4cHh6OtWvX5rmMK+xXTG7ysGTJEgwdOhTjx4/H3r170aRJE0RERCA+Pt7q/Dt37sSzzz6Lfv36Yd++fejZsyd69uyJQ4cOFXHkRc/ebQVIV92XL1/OGs6ePVuEERsnOTkZTZo0wezZs/M1/5kzZ9CjRw906tQJ0dHRGDx4MPr374/169c7OVLj2butNMeOHbPYtypXruykCF3D1q1bMWjQIOzatQuRkZFIT09Ht27dkJycbHOZkny8Ksj2AkrmMSs/JktqAAAFYklEQVQkJATTpk3Dnj17sHv3bnTu3Bn/+te/cPjwYavzu8x+pcim1q1bq0GDBmU9z8zMVNWqVVNTp061Ov/TTz+tevToYTGtTZs2auDAgU6N0xXYu63mz5+vAgMDiyo8lwVArVy5Ms95RowYoRo2bGgxrXfv3ioiIsKZobmc/GyrzZs3KwDq+vXrRRSVa4qPj1cA1NatW23OU5KPVznlZ3vxmKUrX768+vrrr62+5ir7FUtubEhLS8OePXvQtWvXrGkeHh7o2rUroqKirC4TFRVlMT8ARERE2JzfXRRkWwHAzZs3ERYWhtDQ0DzPBEq6krpfFUbTpk1RtWpVPPjgg9ixY4fR4RS5xMREAECFChVszsP9Spef7QXwmJWZmYnFixcjOTkZ4eHhVudxlf2KyY0N//zzDzIzM1GlShWL6VWqVLF5/T42Ntau+d1FQbZV3bp1MW/ePKxevRrfffcdzGYz2rZtiwsXLhRFyMWKrf0qKSkJt2/fNigq11S1alXMmTMHK1aswIoVKxAaGoqOHTti7969RodWZMxmMwYPHoz7778fjRo1sjlfST1e5ZTf7VWSj1kHDx5EmTJl4Ovri1dffRUrV65EgwYNrM7rKvuVa9wEgkqc8PBwi8y/bdu2qF+/PubOnYv33nvPwMioOKtbty7q1q2b9bxt27Y4deoUZs6ciYULFxoYWdEZNGgQDh06hN9//93oUIqF/G6vknzMqlu3LqKjo5GYmIjly5ejT58+2Lp1q80ExxWw5MaGSpUqwdPTE3FxcRbT4+LiEBwcbHWZ4OBgu+Z3FwXZVjl5e3ujWbNmOHnypDNCLNZs7VcBAQEoVaqUQVEVH61bty4x+9Ubb7yBX375BZs3b0ZISEie85bU41V29myvnErSMcvHxwe1a9dGixYtMHXqVDRp0gSffvqp1XldZb9icmODj48PWrRogY0bN2ZNM5vN2Lhxo81rjeHh4RbzA0BkZKTN+d1FQbZVTpmZmTh48CCqVq3qrDCLrZK6XzlKdHS02+9XSim88cYbWLlyJTZt2oQaNWrccZmSvF8VZHvlVJKPWWazGampqVZfc5n9qkirLxczixcvVr6+vmrBggXqyJEjasCAAapcuXIqNjZWKaXUCy+8oEaNGpU1/44dO5SXl5f66KOPVExMjBo/frzy9vZWBw8eNOojFBl7t9XEiRPV+vXr1alTp9SePXvUM888o/z8/NThw4eN+ghF5saNG2rfvn1q3759CoD6+OOP1b59+9TZs2eVUkqNGjVKvfDCC1nznz59Wvn7+6vhw4ermJgYNXv2bOXp6anWrVtn1EcoMvZuq5kzZ6pVq1apEydOqIMHD6q3335beXh4qN9++82oj1AkXnvtNRUYGKi2bNmiLl++nDXcunUrax4er3QF2V4l9Zg1atQotXXrVnXmzBl14MABNWrUKGUymdSGDRuUUq67XzG5uYNZs2apu+++W/n4+KjWrVurXbt2Zb3WoUMH1adPH4v5ly5dqu655x7l4+OjGjZsqNasWVPEERvHnm01ePDgrHmrVKmiHn74YbV3714Doi56WnPlnIO2ffr06aM6dOiQa5mmTZsqHx8fVbNmTTV//vwij9sI9m6rDz74QNWqVUv5+fmpChUqqI4dO6pNmzYZE3wRsraNAFjsJzxe6QqyvUrqMevll19WYWFhysfHRwUFBakuXbpkJTZKue5+ZVJKqaIrJyIiIiJyLta5ISIiIrfC5IaIiIjcCpMbIiIicitMboiIiMitMLkhIiIit8LkhoiIiNwKkxsiIiJyK0xuiIiIyK0wuSEiIiK3wuSGiIiI3AqTGyIiInIr/w80PgwBvUhFLwAAAABJRU5ErkJggg==",
      "text/plain": [
       "PyPlot.Figure(PyObject <matplotlib.figure.Figure object at 0x135227210>)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "PyObject <matplotlib.legend.Legend object at 0x145e6ae10>"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "PyPlot.pygui(false)\n",
    "x = collect(1:length(Ep)) * dt\n",
    "plot(x, Ep, \"b-\", x, Ek, \"g-\", x, Ep+Ek, \"r-\")\n",
    "legend((\"E_pot\", \"E_kin\", \"E_tot\"))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Julia 0.6.2-pre",
   "language": "julia",
   "name": "julia-0.6"
  },
  "language_info": {
   "file_extension": ".jl",
   "mimetype": "application/julia",
   "name": "julia",
   "version": "0.6.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}
